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PREFACE 


ilils document contains the proceedings of the Air Force/NASA Workshop on 
Modeling, Analysis, and Optimization Issues for Large Space Structures held in 
Williamsburg, Virginia, May 13-14, 1982 The workshop was Jointly sponsored by 
NASA Langley Research Center, the Air ForcV^^fflce of Scientific Research, and 
the Air Force Wright Aeronautical Laboratories >^The theme of the workshop was 
modeling, analysis, and optimization of large space structures. Including structure- 
control interaction. Speakers were drawn primarily from Industry, with partici- 
pation from universities and government. The workshop was organized Into three 
sessions: mathematical modeling, analysis methodology, and optimization for 

controllability. Results of the workshop were discussed In a final session .' — ^ r /] ‘ 
Sunmarles of the sessions were presented by session technical secretaries, and ^ 
general discussion followed. In addition to this fourth session, ample time was 
allowed within each session for discussions on topics of Individual papers. Intro- 
ductory remarks were made by Dr. Michael J. Selklnd, Air Force Office of Scientific 
Research, and Mr. Robert C. Goetz, NASA Langley Research Center. A list of workshop 
attendees Is Included In the front of this document. 


The workshop organizers express their appreciation to the session chairmen, 
speakers, and panelists, whose efforts contributed to the technical excellence of 
the workshop. Session chairmen were Dr. Larry D. Pinson, NASA Langley Research 
Center; Mr. M. A. Ostgaard, Air Force Wright Aeronautical Laboratories, and Dr. 

V. B. Venkayya, Air Force Wright Aeronautical Laboratories. Thanks are also due 
to session secretaries Dr. J. Housner, Dr. John Cubser, Dr. V. B. Venkayya, and 
Mr. B. Hanks. 


Larry D. Pinson 

NASA Langley Research Center 

Anthony K. Amos 

Air Force Office of Scientific Research 
V. B. Venkayya 

Air Force Wright Aeronautical Laboratories 
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INTRODUCTION 


V 



^ The workshop presentations ranged over many topics in large space structures, 
including structure-control interaction, structural and structural dynamics 
modeling, ■ thermal analysis, testing, design, and optimisation. The interdisciplinary 
area of structure-control interaction, which is a challenge to analysts, designers, 
and test engineers, was clearly emphasized. Not addre: sed in the workshop was the 
important subject of structural deployment. 

Structure-control interaction is emerging as a separate discipline in which 
structural dynamicists and controls engineers each become proficient in both areas 
to accomplish viable designs. Performance requirements dictate an integrated 
approach to design. The necessity for this merging of structural dynamics and 
controls disciplines has become apparent because of extreme requirements in potential 
large space structures applications. Furthermore, research is now ongoing in 
aeronautics in the areas of active flutter suppression, load alleviation, and 
reduced static stability, all of which require consideration of the control of 
structural motion at a rate that excites structural vibration. Presentations on 
this subject at the workshop revealed the immaturity of the technology. Theoretical 
considerations dominate research. Most controllers are considered to be ideal, and 
the effects of various kinds of structural or actuator nonlinearity are not known. 

The presentation by Lyons and Aubrun of results from experiments on elementary 
structures with a few sensor— actuator combinations showed some agreement between 
analysis and test data. Although other experimental work is known to be under way, 
a reasonable assessment is that several theoretical approaches exist in this area 
with insufficient experience in their application to form a blending of these 
methods which can be applied confidently in a practical situation. Only through 
ground test and analysis programs Involving relatively complex structures and 
associated control systems and through careful space flight experiments will this 
necessary confidence and experience be achieved. Such programs also will cause 
appropriate organizational realignments to enhance communications, as well as 
appropriate merging of the structural dynamics and controls disciplines. 


The trends in research in structural analysis and modeling are toward 
increasing finesse to achieve efficiency. These trends are driven by two consider- 
ations: (1) potential large space structures have too many elements to model with 

the conventional finite element approach, and (2) small, accurate models are a 
necessity for practical design procedures involving resizing and repeated analyses. 
One trend is toward models in vhich structures with many repeating elements are 
idealized as equivalent contlnua, followed by analyses that are a blend of clr.ssical 
partial differential equation solutions and approximate techniques. Assumed- 
function approaches that extend classical Rltz-type methods to nonlinear responses 
of these idealized structures are being applied with success. None of these new 
techniques, however, has been incorporated into generally available computer 
analysis programs. In addition, experience with new structural concepts such as 
cable-stiffened structures is sparse. With extended experience using large space 
structure concepts and Incorporation of new techniques into practical computer 
analysis programs, structural analy-^ls and modeling of most forseeable large space 
structures will become feasible. 

Ground testing is perceived as a major challenge In large space structures 
applications because of gravitational and atmospheric effects. To derive maximum 
benefit from ground tests, it is necessary to have flight data to substantiate 
analytical corrections for suspension system, gravitational, and atmospheric effects. 


Controversy exists concerning the nature of such flight data. Research-oriented 
engineers tend to favor basic, relatively Inexpensive experiments with specific 
objectives, whereas project and systems engineers opt for more complex specimens 
with much greater potential for practical demonstration. As was revealed at the 
workshop, ground test data from tests of large structures configured for space 
applications are scarce, and no data exist for space tests. Such data are needed 
to obtain confidence In the analyses necessary to justify commitment to real 
applications. 

The emphasis In thermal analysis is on obtaining greater analysis efficiency. 
Very little is being done to verify analytical approaches or to ascertain sensi- 
tivities through correlation with tests. Novel approaches to problems such as 
Interelement radiation were revealed. These approaches involve consideration of th 
probable importance of various geometric effects. Assumed-function approaches are 
also being developed to enhance efficiency. 

Examples of the application of structural optimization procedures In the 
design of complex structures are rare. Subassemblies can be designed to static 
load requirements, but the technology for design for dynamic loading Is immature. 
The problem Is complex, and research has been under way for over a decade. Dis- 
cussions at the workshop Indicated that real applications will not occur until 
mission requirements force such a formal approach to design. Large space structuri 
applications seem to provide these requirements. 
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SESSION 1 



MATHEMATICAL MODELING 

Chairman: L. D. Pinson 

Secretary; J. M. Houc-ner 
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ASSESSMENT OF CURRENT STATE OF THE ART IN MODELING 
TECHNIQUES AND ANALYSIS METHODS FOR LARGE SPACE STRUCTURES 


Ahmed K. Noor 

George Washington University Center 
NASA Langley Research Center 
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INTRODUCTION 


Needless to say, a critical and In-depth review of the state of the art In 
modeling techniques and analysis methods for large space structures would require 
much more space than can be allotted to this paper. Therefore, this paper focuses 
on certain aspects of the subject pertaining to the structures discipline. Other 
disciplines such as thermal analysis and modeling and controls are not covered 
In this paper. 

Figure 1 shows examples of the five categories of large space structures sug- 
gested for various applications (see, for example. Refs. 1, 2 and 3). They Include 

1) Booms and other one-dlmenslonal configurations 

2) Planar surfaces 

3) Antennas and curved structures 

A) Platforms 

S) Space stations 

Because of their low mass and high stiffness, repetitive lattice trusses have been 
selected as the primary candidates (Refs. 4 and 5) for most of the large space 
structures. This paper focuses on modeling techniques and analysis methods for 
these structures. 





SPACE nATWXM 

(a) 

Figure 1 









CURRENTLY USED APPROACKES FOR ANALYZING 
REPETITIVE LATTICES 

A review of the state of the art In the analysis, design, and construction of 
lattice structures until 1976 is given In Refs. 6 and ?. The currently used ap- 
proaches for analyzing large repetitive lattices can be grouped into four classes 
as shown in Figure 2. 

The first approach Is the direct method wherein the structure is analyzed as a 
system of discrete finite elements. It has the obvious drawback of being computa- 
tionally expensive for large lattices. TTiis is particularly true when a buckling, 
vibration, or a nonlinear analysis Is required. To remedy this drawback, techniques 
were developed for substantially reducing the number of degrees of freedom in the 
buckling and nonlinear analyses (Refs. 8 and 9). 

The second approach Is based on replacing the actual lattice by a substitute 
continuum model (see, for example. Refs. 10 through 14). It has the limitations 
that: a) the local deformation effects are typically not accounted for, and b) ordi- 

nary (or classical) continuum is not suitable for lattices with rigid or flexible 
joints. To overcome these drawbacks, continuum models have been developed which 
Include the local deformation modes (Refs. 15, 16 and 17). For lattices with rigid 
Joints, mlcropolar continue have been developed (Refs. 18 through 21). 

The third approach is the discrete field method which takes advantage of the 
regularity of the lattice and is based on writing the equilibrium and compatibility 
equations at a typical joint and then using the Taylor series expansion to replace 
these equations by differential equations (see, for example. Refs. 22 and 23). This 
approach works well for simple lattice configurations, but becomes quite involved for 
lattices with complex geometry. 

r ^ C 

The fourth group of methods Is called the periodic structure approach. 

This approach is based on either: a) the combined use of finite elements and transfer 

matrix methods, which is efficient only for rotationally symmetric or simple 
geometries (Ref. 24), or b) the exact representation of the stiffness of an individual 
member from which the analysis of beam-like lattices with simply supported edges can 
be performed (Refs. 25 and 26). The limitations of the periodic structure approach 
can be removed by combining this approach with the substitute continuum approach. 

This paper will focus on the reduction methods and Improved continuum models. 



CURRENTLY USED APPROACHES FOR ANALYZING LARGE REPETITIVE LATTICES 


APPROACH 

DRAWBACKS 

POSSIBLE REMEDIAL ACTIONS 

• DIRECT METHOD 

• STRUCTURE ANALYZED AS 
A SYSTEM OF DISCRETE 
FINITE ELEMENTS 

• COMPUIATIONALLY VERY 
EXPENSIVE FOR LARGE 
LAHICES 

• USE OF REDUCTION MHUODS 
FOR BUCKLING AND NONLINEAR 
. PROBLEMS 

• SUBSTITUTE CONTINUUM 
APPROACH 

• LOCAL DEFORMATION EFFECTS 
NOT ACCOUNTED FOR 

• ORDINARY CONTINUUM NOT 
SUITABLE FOR LAHICES 
WITH RIGID JOINTS 

♦ INCLUDE LOCAL DEroRMATION 
MODES IN THE DEVELOPMENT 

•' USE Ml CROPOLAR CONTI NUA 

. ■ ‘ T-- * ? ^ ''“t ’ 

• DISCRETE FIELD METHOD 

• EQUILIBRIUM AND COMPATIBILITY 
EQUATIONS AREWRIHENAT 

A TYPICAL JOINT 

• TAYLOR SERIES EXPANSIONS 
ARE USED TO DEVELOP 
DIFFERENTIAL EQUATIONS 

• DEVELOPMENT AND SOLUTION 
CAN BE SUBSTATtlTIAL FOR 
COMPLEX LAHICE 
CONFIGURATIONS 


• PERIODIC STRUCTURE 
APPROACH 

• SUBSTRUCTURING - COM- 
BINED USE OF FINITE 
ELEMENTS AND TRANSFER 
MATRIX METHODS 

• EXACT REPRESENTATION 
OF STIFFNESS OF IN- 
DIVIDUAL MEMBERS 

• NOT EFFICIENT FOR COMPLI- 
CATED CONFIGURATION 
AND/OR FOR TRANSIENT 
RESPONSE CALCULATIONS 

• LIMITED TO BEAMLIKE 
LAHICES WITH SIMPLY 
SUPPORTED ENDS OR RING 
CONFIGURATION WITH RIGID 
CENTRAL MAST 

• COMBINE WITH SUBSTITUTE 
CONTINUUM APPROACH 


Figure 2 



OBJECTIVES AND SCOPE 


The objectives of this paper are listed In Figure 3. They are: 

1> To review recent progress In continuum modeling and reduction methods which 
are applicable to large space structures 

2) To Identify some Of the analysis and modeling needs for future large space 
structures 

This paper Is divided Into three parts. The first part deals with continuum 
modeling. Both beam^llke and plate-llke lattices are considered. Linear thermo- 
elastic static response, free vibrations, and buckling problems are treated. The 
lattices can have either pin or rigid Joints. Continuum models have also been 
developed for beam-llke lattices with open thin-walled section longerons, but will 
not be considered in this paper. The second part of the paper deals with reduction 
methods as applied to bifurcation buckling, nonlinear static, and dynamic responses. 
The third part deals with analysis and modeling needs. 


OBJECTIVES 

• REVIEW RECENT PROGRESS IN CONTINUUM MODELING AND REDUCTION 
METHODS WHICH ARE APPLICABLE TO LARGE SPACE STRUCTURES. 

• IDENTIFY ANALYSIS AND MODELING NEEDS. 

SCOPE 

CONTINUUM MODELING 

• BEAM-LIKE AND PLATE-LIKE LAniCES WITH DIFFERENT CONFIGURATIONS. 

• LINEAR THERMOELASTIC STATIC RESPONSE. FREE VIBRATIONS AND 
BUCKLING PROBLEMS. 

• LAniCES WITH PIN AND RIGID JOINTS. 

• BEAM-LIKE LAniCES WITH OPEN THIN-WALLED SECTION LONGERONS. 
REDUCTION METHOD S 

• BIFURCATION BUCKLING, NONLINEAR STATIC AND DYNAMIC RESPONSES. 

ANALYSIS AND MODELING NEEDS 

• LOADS DETERMINATION 

• MODELING AND NONCLASSICAL BEHAVIOR CHARACTERISTICS. 

• COMPUTATIONAL MODELS AND ALGORITHMS. 

Figure 3 


EFFECTIVE THERMOELASTIC CONTINUUM MODEL 

The characteristics of an effective continuum model are outlined In Figure 4 for 
a typical double-layered grid such as the one shown In the figure. It Is a continuum 
having the same amounts of strain and kinetic energies as the original lattice when 
both are deformed Identically. The temperature distribution, loading, and boundary 
conditions simulate those of the original double-layered lattice grid. 

The original three-dimensional double-layered lattice is replaced by a two- 
dimensional continuum plate model. The last two characteristics are perhaps the most 
Important In terms of new developments. These characteristics are 

1) Local deformations are accounted for 

2) Lattices with rigid jolucs are modeled by mlcropolar contlnua 


• ISA CONTINUUM WHICH HAS SAME AMOUNT OFTHERMOEIASTIC STRAIN 
ENERGY STORED IN IT AS ORIGINAL DOUBLE-LAYERED GRID WHEN BOTH 
ARE DEFORMED IDENTICALLY. 

• TEMPERATURE DISTRIBUTION AND BOUNDARY CONDITIONS OF CONTINUUM 
SIMULATE THOSE OF DOUBLE-LAYERED LATTICE GRID 

• THREE-DIMENSIONAL DOUBLE-LAYERED lAHICE IS REPLACED BY TWO- 
DIMENSIONAL CONTINUUM PLATE MODEL. 

• LOCAL DEFORMATIONS "ACCOUNTED FOR 


• ORDINARY CONTINUUM MODELS FOR LAHICES WITH PIN JOINTS AND 
MICROPOLAR CONTINUUM MODELS FOR LATTICES WITH RIGID JOINTS 



LATTICE STRUCTURE TWO-DIMENSIONAL CONTINUUM MODEL 


IjOCAL DEFOBM&TIOHS in AXIM.LY LOADED PLANAR TRUSS 



The local deformation In an axially loaded planar truss is shown in Figure 5. 

The actual deformation has the zig-zag pattern shown on the top sketch. On the 
average, however, the chord members remain straight. The classical continuum averages 
these deformations as shown in the bottom sketch thereby substantially overestimating 
the axial stiffness. The continuxim models developed in Refs. 15, 16 and 17 do account 
for local deformations such as the one shown here. 
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ACTUAL DEFORMATIONS 


DEFORMATION PREDICTED BY CLASSICAL CONTINUUM MODEL 

Figure 5 





ORDINARY VERSUS MICROPOLAR CONTINUA 

A contrast between the ordinary and mlcropolar contlnua Is made In Figure 6. 

For an axially loaded pln-Jolnted truss member the transverse motion is completely 
characterized by the joint displacements. The member rotation is related to the 
joint displacements Wj^ and Wj. Therefore, the appropriate continuum to use in 
modeling pin-jointed trusses is the ordinary continuum for which the dil>p-ta.CCmznt 
£-ceZd comptzteZy cha^cteAizu th& motion o£ the. itfiactuAe. 

On the other hand, for a rigid-jointed member, the transverse motion is charac- 
terized by both the joint displacements w^, as well as the joint rotations 6^, 6j 
which are Independent of the displacements. Therefore, the appropriate continuum to 
use in modeling rigid-jointed flexural members is one whose motion is choAacteAtzed 
by both a dCspiacejnent £teZd and an independent Kozution £ietd. The mlcropolar con- 
tinuum is svtch a continuum (Ref. 27). 


TRUSS MEMBER (PIN JOINTS) BEAM MEMBERIRIGID JOINTS) 


DEFORMATION 

PAHERN 


ROTATIONS 


-p (w. - w.) 
L j I 


^ ~ r * 

0j. 0. (JOINT ROTATIONS) 


APPROPRIATE 


ORDINARY 


CONTINUUM <D«SPIACEMENT FIELD ONLY) 


MICROPOLAR 

(INDEPENDENT DISPLACEMENT 
AND ROTATION FIELDS) 


Figure 6 














BEAM-LIKE LATTICES CONSIDERED IN PRESENT STUDY 

Some typical configurations for the beam-like and plate-like lattice trusses are 
shown in Figure 7. The characteristics of the continuum models for these lattices 
are given in Refs. 15 and 16. 
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(a) TRIANGULAR CROSS SECTION 









(b) RECTANGULAR CROSS SECTION 
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FREE VIBRATIONS OF CANTILEVERED DOUBLE-LACED BEAM 


Consider the free vibrations of the cantilevered beam-like lattice with pin . 
joints shovm in Figure 8. The continuum model for this structure is a cantilever 
beam. The stiffness and the mass characteristics of this beam are given In Ref . 15 . 

The accuracy of the lowest seven frequencies obtained by the continuum model is 
shown in the figure for the two cases of five and twenty bays. The exact frequencies 
were obtained by a direct finite element analysis of the actual structure. For five 
bays the maximum error in the third bending frequency is 7 percent and reduces to 
less than 2 percent for twenty bays. As to be expected, the accuracy of the predictions 
of the continuum model increases with the increase in the number of bays. 


n=5 

FIVE REPEATING ELEMENTS 
bi^ith BENDING MODE 
ti = ith TORSIONAL MODE 



n = 20 

TWENTY REPEATING ELEMENTS 
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FREE VIBRATIONS OF A SIMPLY SUPPORTED HEXAHEDRAL GRID 


The second problem considered is that of- the free vibration response of a double- 
layered hexahedral lattice grid with pin joints shown in Figure 9. The continuum 
model for this grid is taken to be a square plate whose stiffness and dynamic charac- 
teristics are given in Ref. 15. 

In order to amplify the effect of local deformations, the areas of the core 
members of the grid were assumed to be twenty times the areas of the core members. 

The accuracy of the lowest six frequencies is shown in the figure. The solid lines 
refer to the continuum which Includes local deformations, and the hatched lines are 
for the case when local deformations are neglected. The importance of Including the 
local deformation is obvious. 



AREA OF CORE MEMBERS 
AREA OF SURFACE MEMBERS ~ 



Figure 9 



SPATIAL VIERENDEEL GIRDER 

This set of results Includes the free vibrations and buckling of beam-like lattices 
with rigid joints. As an extreme case consider the spatial vlerendeel girder shown in 
Figure 10. The continuum model for this structure is a mlcropolar beam whose stiffness 
and dynamic characteristics are given in Ref. 21. 



Figure 10 


FREE VIBRATIONS OF CANTILEVERED VIERENDEEL GIRDER 


The accuracy of the lowest six vibration frequencies obtained by the mlcropolar 
continuum models is shown In Figure 11. Two cases are considered, namely, five and 
twenty bays. For the case of five bays, the maximum error in the third bending fre- 
quency is less than 5 percent and for twenty bays the error is well within 0.1 percent 


\ 

Accuracy of Mlcropolar Continuum Solution 
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mode shapes for cantilevered vierendeel girder 

between flewerel end tnS;„ref sho^! ?' elt.rnete 
with the exact frequencies obtained by the direct flnite**!?**^*^^^^"* given along 
actual lattice structure (shown betweL pere«Ses)? analysis of the 
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Figure 12 
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BUCKLING OF CANTILEVERED VIERENDEEL GIRDER 







The lowest two bifurcation buckling modes of a ten-bay cantilevered vierendeel 
girder are shown In Figure 13,. The predictions of the micropolar beam model are 
given along with those obtained by direct finite element analysis of the actual 
lattice (shown between parentheses). 


Ten Repeating Elements 




Figure 13 
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REDUCTION METHODS FOR BIFURCATION 
BUCKLING AND NONLINEAR ANALYSIS 


The second topic considered in this paper is the use of reduction methods for 
bifurcation buckling and nonlinear analysis. The basic features of reduction methods 
are outlined in Figure lA. They are techniques for reducing the number of degrees of 
freedom through the transformation shown. The vector {xl represents the original 
displacement degrees of freedom. The vector refers to amplitudes of displacement 
modes and [F] is a transformation matrix whose columns represent a priori chosen global 
displacement modes. 


As to be expected, the effectiveness of reduction methods depends to a great 
extent on the proper selection of displacement modes. In a number of studies, it was 
shown that an effective choice of the displacement modes includes the various order 
derivatives of the displacement vector with respect to the load parameter (Refs. 8 
and 9). These vectors are g&neAated bij (Ulng the. iltuXe etement model of the OKi- 

3X 

gtnal lattice ^tAuctuae. The recursion formulas for evaluating the derivatives 
2 3 ^ 

{-5-^}, ... are obtained by successive differentiation of the original finite 

3p^ 3p^ 

element equations. The left-hand sides of the recursion formulas are the same (see 
Ref. 8). Therefore, onlif one mat/Ux iactofoLzation AequlAed ^OA the geneAation 
o(i all global £unctlx)ni. Several numerical experiments have demonstrated the 
effectiveness of this choice (see Ref. 9). 


DEFINITION: ARE TECHNIQUES FOR REDUCING THE NUMBER OF D.O.F. 
THROUGH THE TRANSFORMATION 





r < < n 


1X1 = original DISPIACMENT D.O.F. IN THE FINITE 
ELEMENT MODEL 

rr] = MATRIX OF GLOBAL DISPLACEMENT MODES 

... = REDUCED D.O.F. - AMPLITUDES OF DISPLACEMENT 
»'* MODES 

SELECTION OF GLOBAL DISPLACEMENT MOOES 



p = LOAD PARAMETER 


• COLUMNS OFCn GENERATED BY USING THE FINITE ELEMENT MODEL OF THE 
ORIGINAL STRUCTURE 

• THEIR GENERATION REQUIRES ONLYONEMATRIX FACTORIZATION 

• NUMERICAL EXPERIMENTS HAVE DEMONSTRATED THEIR EFFECTIVENESS 

Figure 14 


BASIC EQUATIONS USED IN REDUCTION METHODS FOR BIFURCATION 
BUCKLING AND STATIC NONLINEAR PROBLEMS 

The basic equations used in reduction methods for bifurcation and static non- 
linear problems are given in Figure 15. It is worth noting that the original 
displacement unknowns {X} can be on the order of thousands whereas the reduced 
unknowns {<{(} are typically twenty or less. This is true regardless of the complexity 
of the structure and/or the loading. The details of the computational procedure for 
tracing the load— deflection paths in nonlinear static analysis, including identifi- 
cation of bifurcation and limit points, are given in Refs. 28 and 29. Application 
of reduction methods to nonlinear dynamic problems is discussed in Ref. 9. 



ACTUAL (LARGE) PROBLEM 

REDUCED (SMALL) PROBLEM 

FUNDAMENTAL 

UNKNOWNS 



1X1= INDIVIDUAL - lX|=[r]|«ir} AMPLITUDES | 

DISPLACEMENTS 
• THOUSANDS OF UNKNOWNS 

Ul- UlbfiALt/VltNI 

MODES 

• TWENTY OR LESS 

GOVERNING 

EQUATIONS 

BIFURCATION 

BUCKLING 

[K] + p[Kg]]lX| = 0 

[CK] + p [Kg]] 1^1= 0 

STATIC 

NONLINEAR 

RESPONSE 

[KliXI + IG (X)l- plPl = 0 
~ 1000 EQUATIONS 

CK]|<(^} + IG (<^r)}-plp| = 0 
~ 20 EQUATIONS 

HOW TO TRACE LOAD- 
DEFLECTION PATH 

• REPEATED SOLUTION OF 
LARGE SYSTEMS OF 
SIMULTANEOUS NONLINEAR 
ALGEBRAIC EQUATIONS 

• GENERATION OFCF] 

• MARCHING WITH SMALL 
SYSTEM OF EQUATIONS 

• ERROR SENSING AND 
CONTROL (UPDATING [F] 
WHENEVER NEEDED ) 


[K3 = [r]T[K][n . [Kg] = [nT[K][n . {g(^)| =[r]^iG(X)i 


Figure 15 




APPLICATION OF REDUCTION METHODS TO BIFURCATION 
BUCKLING OF A LATTICE TRUSS 

As a simple demonstration of the effectiveness of reduction methods to bifurcation 
buckling problems, consider the thirty-bay cantilevered beam-like lattice truss with 
pin joints subjected to axial loading. The truss has 372 displacement degrees of 
freedom. The convergence of the lowest two buckling loads with the Increase In the 
number of global functions or reduced degrees of freedom Is shown In Figure 16. With 
six global functions, l.e., a 6x6 eigenvalue problem, the results obtained by the 
reduced system are Identical to those obtained by the full system of equations to 
four significant digits. 


30-BAY CANTILEVERED BEAM-LIKE LAHICE TRUSS SUBJECTED TO AXIAL LOADING 






ANALYSIS AND MODELING NEEDS FOR LARGE SPACE STRUCTURES 


The third part of the paper deals with the analysis and modeling needs of large 
space structures. A list of some of the areas where more work Is needed In the future 
Is given in Figure 17. The first area In which work is needed Is the general deploy- 
ment/erection analysis capability which Includes the elements of both rigid body 
kinematics and flexible body dynamics. 

The second area Is that of analysis and modeling of structures with very slender 
members. This Includes members with length to radius of gyration ratio of the order 
of 1,000. It also Includes tension-stiffened and tension-stabilized structures. 

Considerably more work is needed on sensitivity analysis, that is the sensitivity 
of the response of the structure to the various items listed in Figure 17. 

Then there Is a need for a hybrid continuum/discrete modeling and analysis 
capability for specialized problems such as stress concentrations. 

Finally, accurate determination of thermal, dynamic, and control forces is needed 


• CENERAl OCPlOYMCNT/eRECTION ANALYSIS CAPABILITY 

•COMBINED KINEMATIC-STRUCTURAL ANALYSIS 'INCLUDING ROTATION. 
NONLINEAR AND DAMPING CHARACTERISTICS Of JOINTS, AND 
FLEXIBILITY OF MCMBERSI 

• DETERMINATION OF DEPLOYMENT LOADS (BOTH MECHANICAL AND 
OYNAMICI 

• ANALYSIS AND MODELING OF LAHICE STRUCTURES WITH VERY SLENDER 
MEMBERS AND CABLES 

• DYNAMIC CHARACTERISTICS Of STRUCTURES WITH BUCKLED MEMBERS 

• STUDY OF NONE I WAR (LARGE ROTATION) EFFECTS 

• IMPERFECTION SENSITIVITY 

• LEVEL OF MODELING REQUIRED (TENSION-STABILIZED STRUCTURES 
AND wrinkled MF.MBRAWSI 
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• SENSITIVITY ANALYSIS • SENSITIVITY OF RESPONSE TO; 

• SURFAa INACCURACIES (E.G.. CURVED SURFACE STRUCTURES) 

• VARIATIONS IN DESIGN VARIABLES 'REQUIRED FOR EVALUATION 
OF STRUCTURAL CONCEPTS AND FOR OPTIMIZATION) 

• MODELING DETAILS 

• JOINT REXIBILITIES AND DAAIPINC CHARACTERISTICS 

• FAILURE OF SOME MEMBERS (DAMAGE TOLERANCE) 

• HYBRID CONTINUUM - DISCRETE MODELING AND ANALYSIS CAPABILITY 

• DETERMINATION OF THERMAL. DYNAMIC LOADS AND CONTROL FORCES 

• THERMAL LOADS INCLUDING COMFUTATION OF RADIATION VIEW- 
FAaORS 

• DYNAMIC LOADS DUE TO DOCKING. MANEUVERING AND ASSEMBLY 

Figure 17 



FUTURE DIRECTIONS FOP. MODELING AND ANALYSIS 
OF LARGE SPACE STRUCTURES 

Af- .ar as future directions are concerned, the driving forces are the needs for 
evaluat on of structural concepts, and accurate prediction of strength, stiffness, and 
fatigue life of large space structures (Figure 18) . There are numerous opportunities 
provided by new advances in computer hardware, firmware, software, CAD/CAM systems, 
computational algorithms and materials technology. 


DRIVING FORCE 

• EVALUATE STRUCTURAL CONCEPTS 

• PREDICT STRENGTH. STIFFNESS. FATIGUE LIFE AND 
DAMAGE TOLERANCE 

• MODELING AND ANALYSIS OF LARGE STRUCTURES 
SUBJECT TO HARSH ENVIRONMENT 

• RELIABILITY AND ERROR ESTIMATES OF SOLUTION 

OPPORTUNITIES - PROVIDED BY NEW ADVANC^ES IN- 

• COMPUTER HARDWARE AND FIRMWARE 

• COMPUTER SOnWARE AND CAD/CAM SYSTEMS 

• COMPUTATIONAL ALGORITHMS 

• MATERIALS TECHNOLOGY 


Figure 18 



RECEKT AND PROJECTED ADVANCES IN COMPUTER 
HARDWARE AND SOFTWARE 

Some of the recent advances In computer hardware and software are listed In 
Figure 19, The new computing systems include the supercomputers such as the CRAY IS 
and the CDC CYBER 205. The opportunities provided by these large computers are 
discussed In Ref. 30. There are minicomputers with the new array processors such as 
the FPS-164 which has 64 bit-length word and a computational speed of the order of 
12 MFLOPS (millions of floating point operations per second). Microprocessors are 
likely to impact large space structures. Then there are the new multiple CPU computers 
like the APPOLLO computers which have several processors running in an asychronous 
manner using the same data base. 

As far as the user- inter face type of hardware is concerned, it is worth mention- 
ing the new verbal (audio) and visual Interfaces. 

Considerable progress has been made in the software area. Software progress 
Includes relational data bases, which can handle large volumes of data and which 
are a product of the IPAD technology at NASA Langley (Ref. 31), and the new geometric 
modeling and graphics systems (such as the AD-2000 and the ANVIL-4000) . Perhaps one 
area which needs more attention is that of artificial intelligence and its exploitation 
in the design of large space structures. 

• NEW COMPUTING SYSTEMS 

• SUPERCOMPUTERS (E.G.. CRAY. CYBER 205) 

• MIN I COMPUTERS/ ARRAY PROCESSORS (FPS-164) 

• MICROPROCESSORS 

• MULTIPLE CPU COMPUTERS IE.G., APPOLLO COMPUTER) - 
SEVERAL PROCESSORS RUNNING IN AN ASYNCHRONOUS 
MANNER USING THE SAME DATA BASE. 

• HARDWARE - USER INTERFACE 

• GRAPHIC DISPLAY 

• VERBAL (AUDIO) AND VISUAL INTERFACES 

• SOFTWARE 

• RELATIONAL DATA BASES (RIM - NASA LANGLEY) 

• GEOMETRIC MODELING AND GRAPHICS SYSTEMS (E.G., AD-2000. 

ANVIL 4000) 

• ARTIFICIAL INTELLIGENCE 

■Figure 19 




IKTEGRATED ANALYSIS AND CAD SYSTEMS 

The efficient design of large space structures requires a strong Interaction 
between a number of disciplines Including structures, controls, and thermal analysis, 
among others. In response to this need. Integrated analysis and CAD systems have 
been developed and are currently upgraded to perform this task. Three examples are 
listed In Figure 20. 


• INTEGRATED THERMAL-STRUCTURAL-CONTROL CAPABILITY 
(BOEING - NASA GODDARD) 

• COMPUTER-AIDED DESIGN SYSTEM 

(GENERAL DYNAMICS, MARTIN MARIEHA - NASA LANGLEY) 

•FUTURE DISTRIBUTED/ INTEGRATED ANALYSIS DESIGN SYSTEM 
AT LANGLEY 


Figure 20 





FUTURE DISTRIBUTED/INTEGRATED ANALYSIS DES1G!^ 

SYSTEM AT LANGLEY 

A schematic of a future distributed/integrated analysis design system at Langley 
Is shovm in Figure *’1. It has an executive and a number of functional modules each 
representing one of the disciplines. These models interact through a common rela- 
tional information manager (RIM) . The planned system, which is already well under 
development, would allow the different modules to be executed on different computers 
(distributed computing). The computers may even be located at different geographic 
locations. 



RELATIONAL 

INFORMATION 

MANAGER 

(RIM) 
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PRIME 
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Figure 21 







SUroiARY 


In summary, three topics are covered in this paper, namely, recent advances in 
continuum modeling, progress in reduction methods, and analysis and modeling needs 
for large space structures (Figure 22) . 

As far as continuum modeling is concerned, an effective and verified analysis 
capability exists for linear thermoelastic stress, birfurcation buckling, and free 
vibration problems of repetitive lattices. However, application of continuum modeling 
to nonlinear analysis needs more development. 

Reduction methods have proven to be very effective for bifurcation buckling and 
static (steady-state) nonlinear analysis. However, more work is needed to realize 
their full potential for nonlinear dynamic and time-dependent problems. 

As far as a'.talysls and modeling needs are concerned, three areas have been 
identified. 

As to be expected, the modeling and analysis of large space structures will be 
strongly impacted by new advances in computer hardware, software, integrated analysis, 
CAD/CAM systems, and materials technology. 


• CONTINUUM MODELING 

• VERIFIED FOR LINEAR STRESS ANALYSIS. BIFURCATION BUCKLING. 
AND FREE VIBRATION PROBLEMS OF REPETITIVE LATTICES 

• NEEDS DEVELOPMENT FOR NONLINEAR ANALYSIS 
® REDUCTION METHODS 

• VERIFIED FOR BIFURCATION BUCKLING AND STATIC NONLINEAR 
PROBLEMS 

• FURTHER DEVELOPMENTS NEEDED FOR NONLINEAR DYNAMIC PROBLEMS 

• ANALYSIS AND MODELING NEEDS 

• GENERAL DEPLOYMENT/ERECTION ANALYSIS 

• LAHICE STRUCTURES WITH VERY SLENDER MEMBERS AND CABLES 

• SENSITIVITY ANALYSIS 

• IMPACT OF NEW ADVANCES IN COMPUTER HARDWARE. SOFWARE AND 
MATERIALS TECHNOLOGY 


Figure 22 
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LATTICE-TYPE SPACE STRUCTURES CAN EXPERIENCE 
ACCUMULATED PENUMBRAL SHADOWING 


Thermal analysis of relatively sparse structures In the space environment has 
customarily omitted consideration of shadowing by up-sun structural members. This 
convention has been frequently questioned In the case of lattice-type structures 
supporting very large, near-planar. Earth-facing surfaces (e.g., antennas). For 
these, significant shadowing can occur whenever the solar vector Is nearly tangent to 
the orbital path. It thus becomes advisable to quantify the shadowing effect, but 
sparse structures present an exceptional element of complexity. A typical sparse 
structural assembly, a parabolic expandable truss antenna (PETA), Is portrayed in 
figure 1. In such assemblies, multiple up-sun (l.e., shadowing) members may yield 
only partial shadowing of an elemental area of interest, with the degree and duration 
of shadowing being a strong function of the size and density of the structural 
assembly. 



Figure 1 


THE TRADITIONAL APPROACH: THERMAL MODELING OF THE ON-ORBIT ASSEMBLY RADIAL ARM 

Space heating analysis procedures usually assume either total shadowing or total 
solar irradiation of individual elemental areas of mapped surfaces. The Vector Sweep 
program (ref. 1) typifies this approach. While the latter convention is reasonably 
sound for structural assemblies of which most -or all geometric dimensions are of the 
same order of magnitude, it is less realistic for sparse structures composed of slen- 
der members having extremely high length-to-diameter (K./d) ratios. Traditional 
methods were employed in modeling an on-orbit assembly (OOA) spacecraft structural 
assembly, as described in reference 2. The analysis model is shown in figure 2. 
Following computation of radiation view factors for spatially oriented arrays of geo- 
metric subelements, incident space-environment heat rates were computed for each sub- 
element at successive points in time throughout the orbit. The heat rate histories 
were then incorporated in a thermal model of the radial-arm subassembly, which yielded 
temperature histories for the radial-arm structural elements. The thermal analyzer 
model nodal arrangement was geometrically Identical to the Vector Sweep program sub- 
element arrangement permitting direct incorporation of the Vector Sweep output heat 
flux histories in the thermal analyzer model. Notwithstanding its level of detail, 
the thermal model is seen to be limited in scope, in that it can address only a local 
region of the total structural assembly. 


• Numbers represent 



Figure 2 


THERMAL RESPONSE AS PREDICTED BY TRADITIONAL METHODS 

On-orbit assembly spacecraft radial-arm cross-member temperature histories 
obtained in the study in reference 2 are shown in figure 3. It will be instructive 
to make a qualitative comparison of figure 3 with temperature prediction data from 
the space structure heating (SSQ) procedure, to be presented later. It is seen 
that only six shadowing members are Implied in the figure 3 temperature histories. 
Thus, accumulated shadowing by multiple, more distant members is not taken into 
account. This is a commentary on the limited geometric sco,/e of the model of 
figure 2. Conversely, the total shadowing (rather than penumbral) inferred 
during the six non-Earth- shadowed intervals is probably excessive. 
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SHADOWING-MEMBER IDENTIFICATION 

The SSQ program permits definition of numbered structural node locations In the 
spacecraft coordinate system, followed by definition of line segments (structural 
members) In terms of bounding node numbers. Slender structural members are assumed 
to be cylindrical, and are assigned values of diameter, thermal mass per unit of 
length, solar absorptance, emlttance, and view factor to space. Partially or totally 
opaque bodies are also defined In terms of bounding line segments, and the latter 
surfaces are assigned transmissivities that may vary as functions of solar angle of 
incidence. 

Following selection of a structural element of interest 1, all spacecraft coor- 
dinate data are transformed (without rotation) to an 1-centered coordinate system. 

The latter data thus define potential shadowing members j. As shown In figure 5, 
the locus of the spacecraft- Sun vector progressing about element 1 through a complete 
orbit describes a cone with Its apex at the 1-centered system origin. Potential 
shadowers are represented as line segments. Simultaneous solution of line and cone 
equations yields sets of intercepts (x, y, z)j which are then examined for residence 
between limiting values (x, y, -)l,j and (x, y, z)2,j. Qualifying sets define Indi- 
vidual members that actually do shadow element of Interest 1, permitting definition 
of the orbital position 0 at which shadowing occurs. At this point, shadowing mem- 
bers j are retained in a new k array. In order of increasing 01 ^. 
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SHADOW ORIENTATION 

It now becomes necessary to examine the orientation of shadows 6^ as viewed 
rom the element of Interest 1. As seen In figure 6, shadows traversing the visible 
olar disc can be classified as being essentially horizontal or nonhorizontal. The 
atter distinction becomes an important discriminator in computing Individual shadow- 
ng durations. The two shadow classifications are compared in figure 6, in which 
ssumed values fc • the radius of and distance to the Sun are shown: 375,735 n.mi. , and 
0,884,432 n.mi., respectively. Of major significance, the end points of a so-cailed 
orlzontal shadow are assumed to enter and leave the visible solar disc, but it Is 
ssumed that the end points of a nonhorizontal shadow are never visible. 


HORIZONTAL 

shadow 


NONHORIZONTAL 

SHADOW 



INC1.INATI0N ANGLE DETERMINES SHADOW CLASSIFICATION 

Development of the test for ^classifying shadow orientation is shown in fig-_ 
ure 7. Unit vector ^ normal to S and Lj^ is obtained by normalizing the product Lj^ x 
S. Incidence angle is the angle described by nj^ and the projection- of on a 
plane normal to 3. Angle 4>i( is employed later in computing the duration of non- 
horizontal shadow Intervals. The scalar projection of on is (L * ft)ki ahd 
(L • n tan divided bv (L • n}^ is the slope of the Lj^ projection normal^to S. As 
shown in figure 7, (C • ^ tan ({>)k divided by the absolute value of vector (Hi + C • n) 
is the tangent of angle This angle is an indicator_of horizontal versus non- 

horizontal shadow orientation; i.e., for = 2 tan“l (Rs/ls|), member k would barely 
qualify as nonhorizontal (see figure 6). At this point, a reasonable assumption is 
Imposed to define a horizontal shadow in terms of Rg/Js|, as shown in figure 7. 
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SHADOWING INTERVAL DURATION 


It now becomes necessary to compute the durations required for each shadow to 
traverse the visible solar disc. Methods employed to compute shadowing Interval dur- 
ations are Illustrated in figure 8, in which the following terms are Introduced: 

” Beginning of Interval 
0E “ End of interval 
A 8 = 0g “85 

d Shadowing member width 
A *• Distance to Intercept at 8 

Horizontal shadow maximum durations are evaluated in terms of the projection of an 
extremely long shadower on a celestial unit sphere. Nonhorizontal shadow durations 
are related to their Inclination angles <{> to the orbit plane, the lower limit being 
twice the bracketed term (a vertlcle shadow). However, the value that Sg and 0g for 
both horizontal and nonhorizontal shadows cannot exceed is always limited by the 
shadower end coordinates. 
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MERCED SHADOWS 


At this point, the k array of shadowing members is ordered In terms of increasing 
0^. Many adjacent shadows (or portions of horizontal shadows) probably merge to form 
a single larger shadow. These conditions must be recfognlzed and adju::ted. The method 
employed for merged horizontal shadows is illustrated in figure 9. Merged nonhori- 
zontal shadows are accommodated by similar procedures. The latter adjustments permit 
creation of a new "m" array of discrete shadowers as would be viewed from the "i" 
element-of-interest , and each shadower is defined in terms of an effective distance 
H, an effective width d, and bounding angles Bg and Og. 





Ll?r a = MAX OF fun'' V + 

V I -i'k.klii 


LET l» = MIN OF sin"' 


•CREATE NEW MEMBi-K kk 

*kk“*k Mkk = “"’‘“‘Mk.k... 

"i:)kk = ‘^^“'"'0k.k.n 


Figure 9 


SIMULTANEOUS NONMERGED SHADOWS 

The SSQ program now creates a new m array of merged and nonraerged shadowing 
members, ordered in terms of increasing 6g. Overlapping intervals represent sub- 
intervals in which more than one distinct shadow is contained within the visible solar 
disc. It is therefore necessary .to create an n array of subintervals of which &0n = 
(0E - QB)n« which the contributing members m are identified. A graphic repre- 

sentation of this procedure is contained in figure 10. 
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INTEGRATION OF SH/'DOW PROJECTIONS ON THE SOLAR DISC 

Penumbral solar heat flux Is computed by integration of projected shadows m on 
the visible solar disc, as portrayed in figure 11. A solar disc area reduction fac- 
tor SHAOqi is computed as shown in figure 11. Horizontal shadows arc currently assumed 
to reside in the center of the solar disc (constant bp,) . For nonhorizontal shadows, 
SHAD,g is computed at eight equally spaced 9 points within each A0„, Interval. Earth- 
shadowing simulation is superimposed by additional solar flux attenuation during 
penumbra passage (assumed linear with 9) and by total eclipse during umbra passage. 
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COMPUTATION OF SHADOWED INCIDENT SOLAR FLUX 

If selected lines j are Identified as bounding an opaque or semiopaque surface, 
only two such lines can be intercepts. During the resulting interval + n ~ 
solar flux Is further attenuated by a constant or variable transmissivity T, or can be 
eclipsed. Solar angle of Incidence Y to a transmissive surface is computed as shown 
in figure 12, in which the two intercepted bounding lines are represented as vectors 
Lj and Cj + and a vector triple product is formed to compute the solar incidence 
angle. Also shown in figure 12 is the total solar flux Incident to element i, but 
attenuated by (1) the summation of SHAD terms, (2) the transmissivity of a semiopaque 
Intervening surface, and (3) angle of incidence to structural element 1, in which S 
is the solar heat flux constant, an input tern. 
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HEAT FLUX SUMMATION AND RESULTING THERMAL RESPONSE 

Incldei t Earth thermal and albedo radiation fluxes are computed in the SSQ pro- 
gram, but shadowing of these diffuse fluxes is not considered. Earth thermal heat 
flux qg is calculated as shown in figure 13, in which p is the Earth albedo factor 
(input), and Fg is a geometric factor for thermal radiation to a cylinder. It is 
contained iii the program as tabular functions of altitude (an input term) and struc- 
tural member angle of inclination to an Earth radius vector. Earth albedo radiation 
is also computed as noted in figure 13, in which 

H ” altitude 

R “ Earth radius 

cos 0g « cos 3 cos 9 

A/rays of solar. Earth thermal, and Earth albedo heat flux incident to structural 
member i are contained in printed output and are also available for graphic display. 

Total absorbed heat flux can be employed in computing the thermal response of an 
elemental length of structural Member of interest 1. A differential energy balance 
is used, in which me is the elemental thermal mass, t is time (a function of 6), q is 
total absorbed heat flux, a is the Stefan-Boltzmann constant, Fg is the element view 
factor for radiation to space, c is the emittance of the element, and a- is its solar 
absorptance. The closed-form solution in temperature T Is evaluated at successive 
points of time. 
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SSQ PROGRAM OUTPUT OPTIONS INCLUDE 
HEAT FLUX AND TEIIPERATURE PLOTS 

The SSQ program is coded la FORTRAN V for the CDC CYBER 172 computer, and can be 
operated In either the batch or Interactive mode. The input format is orderly and 
uncomplicated. Printed output can be limited to heat flux and temperature prediction 
data or can include expanded data identifying shadowing members, merged shadowers, 
and shadowing members (m array) contained in multiple shadowers (n array). DISSPLA 
software is employed in the graphic output subroutine. Tektronix terminal visual and/ 
or hard copy and FR80 hard-copy output are available. Typical FR80 graphic output 
features are shown in figure 14. SSQ program heat flux output can easily be formatted 
as required for input to other thermal- analysis programs. Temperature output can also 
be formatted to accommodate the input requirements of structural-analysis programs 
such as LASS (ref, 4) or NASTRAN. 
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SELECTION OF GENERIC MEMBERS 

It Is possible to achieve substantial simplification of the SSQ shadowing-nember- 
identif Icatlon routine, and we can greatly reduce the task of computing Incident heat 
flux histories. For example, shadowed heat flux Incident to each of 978 members of a 
12-bay PETA structure throughout a complete orbit can be constructed from results com- 
puted for only 108 members. As shown in figure 15, the PETA structure contains nine 
basic member orientations that recur throughout the entire assembly. It is shown in 
figure 5 that the locus of the subvehicle vector about any point throughout a com- 
plete orbit is a cone of half angle (ti/ 2-6), of which 3 is the solar vector angle of 
incidence to the orbital plane. An array of members spanning the entire spacecraft, 
and inclined 3 degrees from the orbital plane, will contain 108 members, 12 each of 
nine generic orientations. Analysis of the members in this array will permit compu- 
tation of shadowed heat flux histories for all members in the entire spacecraft. 
Exceptions will be shadowing effects caused by major external space hardware. As will 
be seen, the latter effects can be computed separately and superimposed on selected 
heat flux histories, where appropriate. 
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CONSTRUCTING A HEAT FLUX HISTORY DATA BASE 

Following identification of a 12-bay transverse array of 108 generic members in 
the PETA structure (fig. 16) and their respective shadowing members., it is proposed 
to compute shadowed incident heat flux histories for five (tentatively) equally 
spaced locations on each of nine generlcally oriented members (fig. 15) of the 12-bay 
generic member array. This data will comprise a heat flux history data file in which 
each heat flux history is identified in terms of (1) one of the nine basic member 
orientations, (2) one of Che five equally spaced member locations, and (3) distance 
from the spacecraft leading and trailing edges. The data file will permit construc- 
tion of heat flux histories for any or all structural members on the spacecraft, as 
illustrated in figure 16. Consider a member residing at dj^, equidistant from the 
leading and trailing edges, as shown in figure 16(b). Incident heat flux to all 
points on this member will have been already assembled in Che data file for the entire 
orbit (0 6 27t). But figure 16(c) portrays a member location, 62 from the leading 

edge and d 3 from Che trailing edge. Incident heat flux Co points on this member will 
be contained in the d2 record for the orbit Interval (0 £ 0 £ ir) , and in the d 3 
record for the orbit Interval (Tf <. 0 £ 2 ti). Heat flux histories for all members can 
thus be constructed rapidly, especially since many members will reside at identical 
locations with respect to the leading and trailing edges of the spacecraft. In Che 
absence of additional major shadowing, such members will have identical heat flux 
(and temperature) histories. 
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MAJOR BODY SHADOWING 

Heat flux attenuation by major shadowing space hardware. If any, must be taken 
Into account before member thermal response can be computed.. Shadows cast by these 
external bodies will traverse the spacecraft surface, and the shadow path is parallel 
to the orbital plane only if B is zero. It will be necessary to Identify all struc- 
tural members any part of which resides within the limits of a shadow path. Major 
shadowing computations can thus be limited to the latter members. A typical shadow- 
ing situation is portrayed in figure 17. Of the shadowing-body coordinate data, the 
maximum and minimum y values correspond to points 2 and 3, respectively. At least 
one end of potentially shadowed members must have a y coordinate no greater than y£. 
For members meeting the latter test, let point 1 be either end of the structural 
member, and 
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If B'max ^ 3 ^ 3*min either end of the member, shadowing of a member will occur. 
Also, if both B' > B on one end of the member while both B* < B on the other end, 
shadowing will exist. Shadowing Intervals for all shadowed points can be computed 
by the methods of figure 8. Major body shadowing effects must now be superimposed 
on the heat flux data records of figure 16. 
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LIMITING THE SCOPE OF TEMPERATURE COMPUTATIONS 

Judicious selection of an orbit for which 0 Is zero (e.g., an equatorial orbit, 

21 March or 21 September) and a preferred PETA attitude can permit temperature com- 
putation for all members of a 12-bay PETA structure in terms of only 163 members 
(plus members subject to major body shadowing). Effects of the choice of PETA atti- 
tude and $ on the task of computing member temperatures are illustrated in figure IS. - 
Dark-shadowed regions of the FETA profiles represent areas in which the temperature 
histories of all members can be generated by computing the temperatures of a single 
transverse array of generic members. Light-shaded areas are those in which each 
member has one, and only one, thermally identical member in the opposite light-shaded 
area. Thus, the PETA orientation of figure 18a (8 ” 0) would require temperature 
history computation for exactly hc?f of the total number of PETA members. The same 
PETA orientation with $ ^ 0 , shown in figure 18b, would yield a limited central region 
In which all temperature histories are represented In a single transverse generic- 
member array, but the above-noted light-shaded temperature correspondence Is limited 
to members situated on the vertical counterllne. All members in the unshaded regions 
of figure 18c would require individual temperature computation. The PETA orientation 
of figure 18d (8 0) /'elds a large central region in which all member temperature 

histories occur in a sl.-gle transverse generic-member array. The latter advantage Is 
largely repeated In the orientation of figure 18e (8 0) . However, the lower light- 

shaded region members of figure 18e will each be thermally identical to those of the 
upper light-shaded region. Temperature histories need only be computed for structural 
members in the shaded region of figure 18e, which will Include an array of generic 
members at the top of the region. Thus, 163 members can thermally model all 978 mem- 
bers in the 12-bay PETA structure. Added to this, as noted earlier, are the 
members that must be subjected to major body shadowing. 
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SSQ IS FAST AND EFFICIENT 

' Assuming a 12-bay FETA analysis application. It Is feasible to create one or 
more data flies containing temperature histories for multiple points on each member 
of the entire PETA structure, although the Individual members will far outnumber the 
temperature histories. The latter records would be' Indexed in terms of (1) types of 
generic members (figure 19), (2) location on the member, and (3) distances from the 
spacecraft leading and trailing edges. Additionally, specific temperature records 
might be Identified solely In terms of specific elements of the structure. It would 
then be possible to create a detailed and accurate transient temperature profile of 
the entire structure by simply constructing a directory assigning identification of 
specific heat flux and temperature records to each member of the structure. To this 
point, a CYBER 720 CP time expenditure of only 3,380 seconds is estimated for the 
entire structure, exclusive of plotting. This data base can be employed as the basis 
for developing rapid-estimation thermal-analysis algorithms, possibly amenable to 
Interactive execution. 


• Analysis of a 6-bay, 240-member PETA required only 
4,14 CP seconds per elemental point (CYBER 720) 

• 108 generic members can model incident heat flux 
histories, & 1 63 members can model temperature 
histories for an entire 1 2-bay, 978-member PETA 
structure. Estimated CYBER 720 CP time is 3,600 seconds 

• Major body shadowing effects require added computation 

• Graphic output requires added expenditure 


Figure 19 




DEVELOPING A RAPID-ESTIMATION ALGORITHM 

It Ic anticipated that examination of the t .ansal-response data base will 
reveal the existence of large groups of structural members that experience 
similar temperature histories. This would permit the analyst to assign the latter 
members to analogous modules for purposes of subsequent simplified analyses. Thus, 
thermal analysis of a single member might reasonably yield the thermal response of 
all members of the module subset. Members subjected to major body shadowing will 
experience incident heat flux histories differing markedly from those of the 108 
generic heat flux schedules of the 12-bay PETA structure, necessita.;lng separate 
Identification and analytical treatment of these members. Therefore, it is probable 
that the location of analogous modules would be much affected by the location of 
major shadows. Moreover, variations in the latter conditions might be expected from 
mission to mission, or at different times In a single mission. Examples of these 
conditions might include the arrival of an upper-stage propulsion system, or the con- 
tinuing reorientation of large solar panels. Thus, If we are working with a specific 
algorithm, and later elect to relocate a large shadowing solar panel, we would prob- 
ably wish to redesign the thermal- analysis algorithm. Application-peculiar shadowing 
effects therefore suggest the possible desirability of partial or total batch execu- 
tion of a baseline thermal analysis as a necessary prerequisite for redesigning ther- 
mal-analyses algorithms. A ‘.’isual search of the graphic data for general temperature 
excursion trends may not be the only means of deducing simplified thermal-analysis 
algorithms. It would be relatively simple to record such statistical data as extreme 
values, mean values, frequencies, and ranges. Whatever the nature of the simplified 
algorithm, the first step In its development must be execution of the baseline SSQ 
code. (See fig. 20.) 


• Identify generic structural members 

• Create an Incident heat flux data base 

• Superimpose major body shadowing effects 

• Create a temperature history data base 

• Create simplified thermal analysis algorithms 
— Examine the temperature data base 

— Identify key members 
— Identify statistical data & trends 


Figure 20 
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LSS CONTROL OBJECTIVES 


The principal objectives for LS£ control are listed below. The general 
objectives range from basic deployment and maneuvering, where some vibration modes 
may be suppressed, to disturbance rejection for very high performance imaging 
applications. The controls selected generally must produce some combination of 
eigenvalue/eigenvector and loads modification in order to achieve the mission 
objectives. 


GENERAL 

• Deployment management and maneuvering 
•Pointing, stabilization, dimensional control 

• Wave front error management (static and d3maraic) 

• Disturbance rejection 

SPECIFIC 

• 0 o o 

• Stiffness (natural frequency) modification 

• Damping augmentation 

e Eigenvector modification 

• Disturbance load modification 

~ Maneuvering 
- Steady state 



STATE-OF-THE-ART ASSESSMENT AND PROBLEMS 

The ultimate mechanization for LSS controls depends on development within at 
least the four areas shown below. The principal difficulties now appear to be a lack 
of funded experimental efforts and inadequate signal processing and actuation hardware. 
The requirements for integrated spacecraft control, including structures, algorithms, 
sensors/actuators, and architectures, are still not generally acknowledged by the 
systems planners. Large high-performance vehicles cannot meet mission performance 
goals without such an approach. 



STUDY AJt£.\ 

1 

1 ACH1EVEHEN7S 

I CURRENT PROBLEMS 

1. 

Deployment end meoeuverlng 

1 

1 - Maneuvering algorithms and 

1 

1 Control requirements and demon- 


control analysis 

* demonstrations 

' stratlons for controlled deploy- 



1 

1 cent largely ignored 



1 - Kultlbody software 

1 



1 - Ocploj-ment ground tests 

1 

2 . 

Syntheels/analysie tools 

1 

, - Basic SAS well understood 

T 

1 - Interface to structural dj-namcs 


and nethods for robust 

1 



control designs 

1 - Robust control design methods 

1 - CAD software for rapid design 



, developed (stability) 




1 

- Robust high performance (low 



1 - Experimental verification 

1 sensitivity) 

3. 

IdentlflcaClon/adaptlve 

— 1 

1 - Algorithms well developed 

i 

1 - Few experimental tests 


Methods for high 




performance 

t - Self-tuning methods applied 

1 - Dependence on FFT method's 



1 to some S/C & helicopter 

1 



* systems 

' - Adaptive methods largely 



1 

1 untested 



1 

- Real-time mechanizations 



1 

1 difficult 

4. 

System archltectnres/ 

1 

1 - A/P systems demonstrated 

1 

1 - No Integrated control philosophy 


mechanizations 

closed-loop (200 Hz) 

* for space structures 

1 



1 - High throughput multiprocessor/ 

1 - Experimental tests still unplanned 



parallel processor architectures 




1 In development 

1 - Existing A/P architectures 



1 

1 Inadequate 



- Hardware tests in progress 




1 

■ 

1 - Inadequate actuator research 

« 



PLATE EXPERIMENT: GENERAL SCHEMATIC 


This picture describes the circular plate suspension and the actuator /sensor 
devices attached to it. The "trapeze" type suspension allows quasi-free motions 
of the plate in sill directions, minimizing the interaction with flexible modes. 
The spring and mass system provides isolation against vibrations coming through 
the mounting bracket. This bracket is attached to a frame bolted on the 14-foot 
optics table. 



PPM ACTUATORS (FRONT SIDE OF PLATE) 
AND ACCELEROMETERS (BACK SIDE) 




ANGULAR OPTICAL SENSING SYSTEM 


This chart shows how the local rotation of the plate is measured at theioirror 
locations. The laser beams (obtained by splitting a single input beam) arc 
reflected by the mirrors and fall upon the 2-axis linear detectors. The (x’ , y’) 
positions of the spots produce pairs of proportional signals. The x' deflection 
showed on the chart corresponds to a rotation of the mirror about the x' axis of 
the plate, and the position y’ corresponds in the same way to a rotation about the y 
axis. 



X 


PLATE 



PLATE EXPERIMENT SENSORS AND ACTUATORS 

This chart depicts the nature and location of the actuators and sensors used 
on the circular plate experiment. The actuators are of two types: 

1. Contactless Actuators — Magnetic forces are exerted on a small magnet 
attached to the plate. This provides the necessary actuation for 
absolute pointing. 

2. Pivoted Proof-Mass Actuators (PPM) - Reaction forces are exerted by 
moving a small mass connected to the plate by an electrodynamlc motor. 

Only AC forces can be produced this way. 

Sensors are of the angular type (deflection of a light beam by a small mirror 
is sensed by a linear photo detector) or of the linear displacement type (measured 
via phase difference of the return light from a corner mirror using microphase 
detectors). In addition, accelerometers are used to measure vibrations along 
the Z axis. 




OVERALL SYSTEM CONFIGURATION 

This chart shows the Interaction and Interfaces between the various processors 
and the test structure. Control gains generated by control synthesis programs on 
the UNIVAC 11—10 are transmitted directly to the STI/DEC 11—23 microprocessor system 
and stored on disk. These gains can then be loaded, when needed, in the Array 
Processor (AP) before starting control experiments. The AP has its own A/Ds and 
DACs and thus carries out the control of the specimen Independently of the rest 
of the system. It can however be directed by the 11-23 to either start or stop 
controlling, or to acquire a time-slice of data for later examination by the 
11-23 software. Also, the 11-23 can concurrently run dynamic characterizations 
of the specimen (either open- or closed— loop) since it can acquire sensor data and 
send excitation signals via its own A/Ds and DACs. 






CONTROL EXPERIMENT PROCEDURE 


This chart is a flowchart shoving the various opera t Iona Involved In the 
structural control experiments. From the finite element analysis, a linear model 
is derived which Is used as a basis for control synthesis programs. These programs 
generate a set of gain matrices which are stored In the 11-23 and may be manually 
updated If needed. The appropriate set Is loaded In the array processor. Then 
the control algorithm is started and the AP controls the specimen and may acquire 
data. This data Is stored on disk and can be analyzed off-line by the 11-23 or 
transmitted to a larger computer (VAX, UNIVAC 11-10) for further processing. Also, 
while the AP Is running, dynamic characterization may be carried out by the 11-23. 










CIRCULAR PLATE RESPONSE TO STEP TORQUE 

This chart shows the angular rotations of the plate about Its x* and y* axes 
when a step torque Is applied about the x' axis. The controller Is digital as 
described In the previous chart. The integral part of the controller makes the 
angle return to zero, exactly cancelling the disturbance (DC) torque. The next 
chart shows the steady-state pointing error. The rms value Is about 5 prad. 


30 uRad 


Transverse axis response 


mi 


ii 
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400 uR«d 


• Direct axis response 


TORQUE OFF | TORQUE ON 



CIRCULAR PLATE: CLOSED-LOOP TEST (HAG 3RB) 

MEASURED AND ESTIMATED OPTICAL SENSOR OUTPUT 

This chart shows comparisons between estimated and measured sensor signals. 
The sensor estimates are obtained by reconstruction from the state estimates. 
Errors In the translational mode estimate are more obvious In the velocity 
comparisons. 
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CIRCULAR PLATE: CLOSED-LOOP TEST (3RB, 2F) 

In this chart, a direct comparison between closed-loop and open-loop responses 1*' 

Is made showing the faster decay Introduced by th.e control system. Also the ^ 

steady-state amplitude has been reduced by a factor of 2. The controller controls 
two flexible modes and uses the two top actuators only. 



INVESTIGATIVE TECHNIQUES 

The least-squares method gives biased estimates of parameters when the input 
measurewnt is contaminated by noise. The instrumental-variables technique was 
developed to mi nimi ze the bias in estimates. In addition, unlike least squares, 
the instrumental-variables method does not require solution to a nonlinear 
programming problem. The instrumental-variable estimates are, therefore, obtained 
in a single step. 

The recursive prediction error methods based on discrete ARHA models update 
estimates at each data point. The algorithms include recursive maximum likelihood, 
recursive least squares, recursive instrumental variables, and extended Kalman filter 
(corrected). These have been used successfully in signal processing and process 
control and monitoring applications. 


ALGQRITHn 


ADVANTAGES 


DISADVANTAGES 


Frequency 

Domoin 

Instrumental 

Varlobles 


•Always converges 

•Requires no starting 
values 

•Data can be averaged 
for use with low 
slgnal-to-noise ratio 

•Number of active modes 
can be easily determined 


•Not minimum variance 

•Can only be block 
recursive 

•Requires solution to 
a polynomial equation 

•Difficult to make It 
robust 


f , 


Recursive 

Masclmuro 

Likelihood 


•Minimum variance 

•Always converges 

•Easy to Incorporate 
robustness against 
data dropouts 


•Requires solution to a 
polynomial equation to 
compute damp* no ratios 
and freouencles 

•Discrete ARMA leads to 
poles near the unit 
clrcTe 
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OPEN-LOOP INPUT AND SPECTRUM FOR PLATE IDENTIFICATION 

The plate open- and closed-loop systems are excited » In these te 
square wave driving a CEM actuator. The corresponding FFT shows that 
excitation energy desirable for system ID Is not Ideal for this exclt 
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CLOSED-LOOP DAMPING RATIO VERSUS MODAL FREQUe’CY 

H..e* damping ratio is plotted versus modal frequency for some of the higher 
frequency modes. Variations from about 30 percent to nearly undamped are observed 
All modeled modes and damping ratios are shown. 
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STRUCTURAL MODELLING AND CONTROL DESIGN UNDER INCOMPLETE 
PARAMETER INFORMATION: THE MAXIMUM-ENTROPY APPROACH 



LSS CONTROL METHODOLOGY REQUIREMENTS 

Figure 1 indicates the overall objectives of the work considered here. Besides 
the obvious basic requirement for a LSS control design methodology, there are several 
additional requirements needed to avoid excessive development cost and possible design 
impasses. Actually, these additional requirements are not demanded merely from a 
desire for mathematical elegance but are often essential for the attainment of the 
.basic objective. Ideally, as the right of the figure Indicates, we should like to 
devise a rigorous, one-step design procedure which . renders inadvertent design of non- 
robust controllers unlikely or Impossible and directly addresses the effects of model- 
ling uncertainty In large-order systems. The present paper concerns a new modelling 
approach which succeeds in attaining many of these desirable features. Put most 
briefly, our approach proceeds by adhering to the following two ground rules: 

1. Incorporate parameter uncertainties directly Into the design model (via 
minimum data/maximum entropy models). 

2. With a hlgh-order stochastic design model, design Implcmentable 
(fixed-order) compensation which is optimal "on the average." (Our 
work has concentrated on mean-square optimization because of Its 
simplicity and relative familiarity). 


AVOIDANCE OF DESIGN IMPASSE 


BASIC REQUIREMENT 

ACCEPTABLE PERFORMANCE 
AT REASONABLE ACTUATION 
LEVELS 



• AD HOC. CUT-ANO-TRY PROCEDURES MINIMIZED OR 
ELIMINATED (Rigorous. One-Step Design 

Procedure) 

• GUARANTEE ELEMENTARY SYSTEM PROPERTIES 
(Robust Stetility) - WITHOUT SERIOUS 
PERFORMANCE DEGRADATION 

(Inedvertent Design 01 Unstable Controller 
Rendered Unlikely Or Impossible) 

• APPLICABLE TO LARGE ORDER SYSTEMS WITH 
MODELLING INACCURACIES 

(Direct Design Under imperfect System 
Information) 


Figure 1 




CONSEQUENCES OF THE CONVENTIONAL MODELLING PHILOSOPHY 

The underlying motivation and ultimate ramifications of the above ideas can be 
explained by contrasting current design strategies with the present approach (Figs. 
2,3). 

Most schemes so far proposed invoke (at least Implicitly) the existence of a 
very high order and allegedly accurate design model which presumes complete informa- 
tion on the values of all structural parameters. As Fig. 2 Indicates, this initial 
modelling step engenders a number of seemingly impenetrable obstacles. The first and 
most serious obstacle, "large dimension", is connected with the difficulties imposed 
upon optimization computations by the large dimension of the plant model. There are 
basically two ways to circumvent this obstacle. The left-hand path shown in Fig. 2 
proceeds by limiting the amount of system information retained in the model, l.e., by 
model simplification employing any one of a number of techniques. Following this 
model reduction step, an "optimal" control design can be computed by means of a suit- 
able LQG-based algorithm. At this point one has a design but no assurance (arising 
from the underlying mathematics of the design formulation) of stability and/or robust 
performance in the face of inevitable spillover effects and parameter errors. Thus, 
it is usually necessary to "hedge" about a nominal LQG-based design In order . to re- 
cover essential system properties - and this involves lengthy and complicated design 
iterations. 

Another approach to the treatment of the dimensionality problem is suggested by 
the right-hand path In Fig. 2. Here, one seeks to limit the system information incor- 
porated in the control by designing a "simple" control of Inherently energy dissipa- 
tive form. This manages to ^.ndle the problems of spillover and other forms of model- 
ling error but Is limited to relatively low performance. The reason is to be found in 
the inescapable trade-off between robustness and performance; one can be "safe" 
everywhere in parameter space but cannot simultaneously achieve excellent performance 
for any particular parameter values. 

Thus, there is difficulty in securing a control which is both robust and optimal 
(in some well-defined sense). Moreover, sine.: the two generic design approaches 
described here do not rigorously guarantee desirable properties, much checking and 
design iteration is required. Finally, each attempt to circumvent a particular 
obstacle seems to lead to further difficulties (for example the problem of large 
dimensionality appears to necessitate modal truncation, but this step, in turn, raises 
the sp. 2 cter of spillover instabilities, etc.). 

The absence of theoretical simplicity and the plethora of ad hoc design steps 
generally indicate chat something is wrong - some fundamental element of realism has 
been Ignored. The view advanced here is that what is at fault is not to be found in 
the detailed design strategies of Fig. 2 (since these are logically driven by the 
initial premise) but rather arises from the presumption of a completely accurate, 
deterministic structural model. Actually, due to numerous sources of modelling error, 
and in advance of extensive testing and Identification, there exists considerable 
a pr*iori uncertainty in the structural parameters. Moreover, one typically has very 
little data with regard to the deviations of parameters from their nominal values. 
Thus, if one rules out the possibility of (off-line) identification of all parameters 
of a very large order system and attempts to devise a control which satisfactorily 
handles most parameter uncertainties at their a pirCori levels, then realism demands 
that such a control design be predicated upon a stochastic model which Incorporates 
the kind of extremely limited statistical information actually available. 


Such a stochastic structural model is one major development of the work reviewed 
herein. In contrast to the customary deterministic model, the minimum data/maximun 
entropy model directly incorporates the least possibla (in the sense to be described 
later) a prion parameter information. The approach taken here is to adopt this 
model as the basic design model (thus incorporating the effects of parameter uncer- 
tainty at a fundamental level), and design mean-square optimal controls (that i«, 
choose the control law to minimize the average of a quadratic performance index oyer 
the parameter ensemble) . 


DETERMINISTIC STRUCTURAL MODEL 
(Praiutnet Complata Information) 


LIMIT INFO. 
RETAINED IN 
MODEL 


LIMIT INFO. 
INCORPORATED 
BY CONTROLLER 


LARGE 


DIMENSION 


SIMPLIFICATION 

TRUNCATION (Aggregation) 
DECOMPOSITION (Perturbation) 


MODEL SIMPLIFICATION 


"OPTIMAL" DESIGN FOR REDUCED 
MODEL 


POSITIVE REAL DESIGN. 
AND/OR RATE-OUTPUT 
FEEDBACK 


SPILLOVER INSTABILITIES 
AND SIMILAR DIFFICULTIES 


MODELLING UNCERTAINTIES 
(Robuatnaaa) 
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PERFORMANCE REQUIREMENTS 
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RAMIFICATIONS OF MAXIMUM-ENTROPY MODELLING 

Figure 3 sketches the major consequences of this approach when applied to the 
typical case in which all but a modest number of system states are associated with 
high levels of parameter uncertainty. First, the maximum entropy modelling philoso- 
phy throws new light on the dimensionality problem. Although the number of 
"poorly known" states may be very large, their statistical response is extremely 
simple and, in particular, the second-moment response of such states may often be 
computed from closed-form expressions. In other words, the process of averaging over 
the parameter ensemble automatically Induces a kind of "stochastic simplification". 
Since we consider averaged quadratic functionals as performance measures, this feature 
extends also to the optimization problem (as will be illustrated below). Thus, the 
control corresponding to the poorly known portion of the system can be determined in 
advance of burdensome computations and a mean-square optimal design is possible for 
very large order systems. 

Ue note that use of the maximum entropy model as the design model gives rise to 
new forms of optimality conditions, representing generalizations of the familiar 
Riccatl equations of LQG. In formulating the optimality conditions for implementable 
compensation, we have proceeded by climbing a ladder; each successive rung is a 
controller form of increasing complexity and realism. 

Even with implementation constraints imposed, the present approach has the poten- 
tial to treat very high order system models. Thus, since spillover effects ultimately 
arise from model truncation and the consequently artificial division of the full-order 
system, "spillover" becomes a phantom within the present context. 

Finally, the essential design model is as unpresumptive of prior information as 
possible and thus envelops the actual parameter statistics in that it Is "maximally 
chaotic". Moreover, the optimization formulation implicitly Imposes robustness con- 
straints, as will be illustrated shortly. Thus, by virtue of the intrinsic mathe- 
matical properties of the original design model, it is possible to secure ah optimal 
and robust design with little or no design It^eration. 


DESIGN MODCl INCORPORATES 
CfAST POSS/Rir PARAMETER 
INFORMATION 


RESPONSE OP **POORLY KNOWN" 
STATES IS YfAr SIMPLE 


MEAN SnUARE OPTIMAL DESIGN 
IS POSSIBLE FOR FUil OfiDi^ 
SYSTEM 


SPILLOVER DOES NOT EXIST 


MINIMUM INFORMATION MODEL 
ENVELOPES ACTUAL 
UNCERTAINTIES 
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LIMITATIONS OF DETERMINISTIC DESIGN: THE NEED TO 


QUANTIFY PRIOR UNCERTAINTY 

Now let us explore the basic ideas of the maximum-entropy modelling approach, 
outline what is Involved in its implementation, and then illustrate the various prop- 
erties described above. 

First, consider a structure whose parameters nay be found in some “parameter 
spacu" as sketched in Fig. 4. a. Standard LQG-based approaches implicitly assume 
that all system maps are known and, consequently, produce a design which Is optimal 
(with respect to a quadratic performance measure) only for a single point (associated 
with the nominal values.) in parameter space. However, due to actual In-misslon 
changes, mathematical modelling errors arising from truncations implicit in the finite 
element method, etc., all system parameters are not , in fact, known. Thus, to put the 
matter in the most general terms, a structural model can never encompass the "truth" - 
rather a model should be regarded as a mathematical statement of what and how much is 
known. Considered as such, a model must not only specify nominal values but must also 
contain an admission of prior Ignorance regarding parameter deviations from expected 
values. An admission of prior uncertainty can be quantified by supposing the param- 
eters to be distributed according to some probability law (see Fig'. 4.b, where the 
shading indicates the region of significant probability). 


MODAL PARAMETER SPACE 



Figure 4. a 
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KIHIKUM EXPECTED COST DESIGN; INCOMPLETENESS OF PRIOR 
PROBABILITIES AND THE MAXIMUM ENTROPY IDEA 

Thus, one might be tempted, as suggested in Fig. A.b, to assume that all neces- 
sary probability distributions are given and proceed to design a control which is 
optimal "on the average" by minimizing the expected value of a quadratic cost. How- 
ever, the problem is more difficult than this because, in practice, a complete prob- 
ability model is never available from empirical determinations. Rather, we need to 
induce a complete probability model from a highly incomplete set of available data. 

A fundamental logical requirement is that this be done in a manner which avoids 
Inventing data which does not exist! In other words, it is necessary to construct a 
complete probability assignment which is consistent with the data at hand but admits 
the greatest possible prior ignorance with regard to all other data* This is the 
heart of the maximum entropy modelling idea. The appropriate quantitative procedure 
has been given by Jaynes (Refs. 1,2): first define a measure of prior ignorance, l.e., 
the entropy ("entropy" as in information theory not thermodynamics), then determine 
the probability law which maximizes this functional subject to the constraints imposed 
by available data. 
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THE MINIMUM DATA/MAXIMUM ENTROPY DESIGN APPROACH 

Having overcome the difficulties Imposed by Incomplete available data. It 
should be noted that for flexible mechanical systems one may Identify a "minimum 
data set" which is just sufficient to induce aqy well-defined maximum entropy model. 
In other words, all admissible sets of available data must include the minimum set, 
and lack of any element of this minimum set will cause the induced maximum entropy 
model to "blow up" in certain crucial respects. Since, in practice, one is provided 
with little or no prior statistical data., it is not only design conservative, but 
also realistic to acknowledge as available data only the minimum data set. 

Thus, as sketched in Fig. A.c, our stochastic design approach involves three 
main stages (see Refs. 3-5). First, the minimum data set is constructed and appro- 
priate numerical values are assigned; next, a maximum-entropy probability model is 
induced from the minimum data (giving the basic design model), and, finally, a mean- 
square optimal design is determined under the maximum-entropy statistics. This pro- 
cedure gives us a mechanism for incorporating incomplete system information within 
the control design. Moreover, as Fig. 4.c suggests, the maximum-entropy model is 
maximally dispersed in parameter space, and one can guarantee that the resulting 
design will very greatly reduce the probability of severe performance degradation 
in the face of parameter deviations. 







MAXIMUM-ENTROPY APPROACH: SPECIFIC FEATURES 

FOR STRUCTURAL SYSTEMS 

Figure 5 indicates more specific features of the maximum-entropy approach as 
applied to linear mechanical systems. In the state-space equation at the top of the 
figure, X denotes the state vector of the plant and any dynamic compensation that 
might be employed, and the dynamics map Is the sum of a deterministic portion, A, rep- 
resenting the nominal parameter values and a zero-mean random portion, a. Hlth this 
form of representation, the formulation Is presently capable of encompassing most of 
the more serious sources of modelling error In structural systems. In particular, a 
may represent general uncertainties in the stiffness matrix and In the placement and 
alignment of actuators and sensors (see Refs. 6 and 7). 

A suitable measure (which Is particularly familiar In information theory) of 
a priori Ignorance Is given In the figure, where p(x) denotes the probability density 
of the state. Actually, since entropy is only a relative measure of information. Its 
definition is not unique, and the expression given below Is one Of several possibili- 
ties (note that Refs. 3 and 5 defined entropy measures directly upon the parameter 
space). Although each possibility Is particularly suited tr specific purpose, all 
these entropy definitions lead to the same results. In part..-ular, besides specifica- 
tion of nominal parameter values and an enumeration of those scalar structural par- 
ameters whose Independent variation is to be Included In a, the minimum data set con- 
sists of the "uncertainty relaxation times". Each of these time scales Is associated 
with a particular uncertain scalar parameter, and is simply and explicitly defined in 
terms of familiar statistical quantities (see Refs. 3 and 5 for details). Qualita- 
tively speaking, the relaxation times give the time scales over which the effects of 
the associated uncertainties upon second-moment response are manifested, and their 
reciprocals constitute fundamental measures of the magnitudes of parameter deviations. 

The maximum-entropy model Induced by this data set is a rather special form of 
Stratonovlch state-dependent noise (see Ref. 8), and It yields a modified Lyapunov 
equation for determination of the state covariance matrix, Q = E(xx^]. The linear 
operator, D[.], appearing at the bottom of Fig. 5 , Is a mapping of the class of 
symmetric matrices onto itself, and It constitutes the principal modification. Since 
the matrix elements of D[Q] are Inversely proportional to the relaxation times, this 
"stochastic Lyapunov equation" reduces to the ordinary Lyapunov equation In the 
deterministic limit (wherein all relaxation times are infinite). Since Q Is the only 



response quantity directly relevant to mean-square performance, optimality conditions 
are obtained In the form of "stochastic Rlccatl equations" involving analogous modi- 
fications of the standard Rlccatl equation. Thus, these new forms of optimality con- 
ditions contain the familiar results for deterministically parametcred systems as a 
special case. 
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EFFECTS OF PABAMETER UNCERTAINTY: INCOHERENCE 

Although Inaccessible to us In practice, a complete probabilistic description of 
the structure nevertheless exists In principle. The maximum-entropiy model may be 
viewed as an approximation to this complete probability law which extracts only cer- 
tain salient features which are to b.r modelled with fidelity. This naturally leads 
to the question: Nhat features of the actual parameter statistics are preserved In 
the maximum entropy model? (Obviously, the relaxation times are preserved, but what 
Is their physical meaning?) 

To provide an answer, let us consider the effects of parameter uncertainty on 
second-moment response of mechanical systems generally - not necessarily under the 
maximum entropy model. First, for purposes of Illustration. Fig. 6 shows a simple 
system of two (weakly) coupled oscillators in which the stiffnesses k. and k, are 
subject to uncertainty. When the magnitudes of the deviations are significant and 
largely Independent of the detailed probability distributions £f k. and k.. the 
cross-correlation coefficient (defined In the figure, where (.) denotes an average 
over the parameter ensemble) of the two oscillators approaches a value of order T 
after a time period of order T following some initial state. Here T Is. the minlirum 
of the relaxation times associated with k and k,. This tendency of uncertainty to 
suppress cross correlation among distinct^modes nas been termed the "incoherence" 
effect, and the situation depicted In Fig. 6 may be generalized In an obvious manner 
to the case of n(n > 2) oscillators. 



T 4 MIN(T,. ' 2) 

T-|. T2 J RELAXATION TIMES OF K, and K2 


Figure 6 



EFFECTS OF PARAMETE.' UNCERTAINTY: ISOTROPY 

Next consider the sane system as In Fig. 6 but with uncertainty only In. the 
coupling spring stiffness (Fig. 7). Besides the Incoherence effect* averaging over 
the parameter ensemble introduces a net power flow between the uncoupled oscillators 
which tends to equalize their separate mechanical energies. For sizable uncertainty 
levels and over a time of order T following some initial state, '_he system relaxes 
to a special statistical equilibrium wherein the relative energy difference Is propor- 
tlonal to T^. Here again, T^ is the relaxation time (as defined In Ref. 5) associated 
with the random portion of k . Thus, for very large uncertainty (T small), this 
"Isotropy" effect produces equality of energies and, from the point*^of view of. their 
second-moment response, the two oscillators are largely indistinguishable. This sort 
of thermal equilibrium Is established In an analogous manner for a large-order system 
wherein uncertainties exist which tend to couple nominally uncoupled modes. 

The effects Just discussed have been known (separately) for some time (see, for 
example. Refs. 9-11) and are real - not artifacts of the maximum -entropy modelling 
approach. Rather, the magltudes and associated time scales of these effects are 
precisely what the maximum- entropy model preserves and faithfully models . 
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IMPACT OF UNCERTAINTIES ON MEAN-SQUARE OPTIMAL DESIGN 


The great Importance of accurate modelling of Incoherence and Isotropy arises 
from the fact that these effects engender distinct qualitative regimes wherein mean- 
square optimal control, design is dramatically affected. To Illustrate this, consider j 
the mean-square optimal design of full-state feedback regulation Crelevant equations J 
are indicated at the top of Fig. 8) for a hlgh-order stochastic model which incorpo- I 

ra tee the "actual" parameter statistics. Also suppose that modes are arranged In ^ 

order of Increasing uncertainty levels. By virtue of the Incoherence effect, there 
emerges an ‘'Incoherent range" wherein the modes are (to a close approximation) mutu^ '] 
ally uncorrelated and uncorrelated with lower order modes. The remaining (low order) 
m^des retain significant cross-correlation and constitute the "coherent range". In . 1 
addition, for sufficient uncertainty In the structural stiffness, there may exist an 
"Isotrcplc subrange" (where approximate equality of modal energies holds) within the ^ 
Incoherent range. Clearly, since the performance Index can be expressed as a linear 
function of the second moments, the contribution of all the incoherent modes Is of the 
form of a weighted sum of mean squares of the separate states (as Indicated In Fig. 8). 
This reduction to a sum of mean squares has a profound effect on the control design I-'h 
because it implies that any control effort devoted to shifting modal frequencies and/or 
closed-loop eigenvectors Is wasted. In consequence, mean-square optimal design demands 
an inherently robust positive real or rate feedback controller for the incoherent modes. 
Obviously, this dramatic qualitative effect cannot be Ignored. 2'j 

Now the maximum entropy model not only models the above qualitative regimes faith- 
fully , but In addition, the associated optimality condition (the stochastic Rlccatl 
equation) automatically produces an Inherently energy dissipative rate feedback control 
for the highly uncertain. Incoherent modes. Moreover, the incoherent range control 
design is in close conformity with the asymptotic solution of the stochastic Rlccatl 
equation for large uncertainty levels. Since, in many instances, this asymptotic 
solution can be given In closed form , one may regard the mean-square optimal design j 

for the very poorly known states as being known in advance of burdensome computations. 
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NUMERICAL EXAMPLE POR THE REGULATOR PROBLEM 

As a specific example, consider the regulator design for a simply supported 
Bernoulll-Euler beam with one force actuator depicted iii Fig. 9. Here an energy 
weighting Is imposed on the state In the performance Index, and we denote by the 
scalar p the relative control Input weighting. (Fjr complete details regarding this 
example, including closed-loop poles, sensitivity studies, etc,, see Ref. 12). Uncer- 
talntles only In the open^loop frequencies are considered, and frequency deviations 
from the nominal values, u , are taken to be Gaussian. The corresponding relaxation 
times, T. , are seen to be Inversely proportional to the relative standard deviations. 
To provide a simple model which reflects progressive degradation In modelling 
accuracy with increasing modal order, we take the to Increase with nominal fre- 
quency. 


EXAMPLE: SIMPLY SUPPORTED BEAM WITH FORCE ACTUATOR 
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DESIGN RESULTS AND QUALITATIVE INTERPRETATIONS - 
STOCHASTIC SIMPLIFICATION 

For specific values of damping, p, and uncertainty level. Fig., 10 shows the mean- 
square optimal design for a 100-mode model. For clarity, we show the gain magnitudes 
scaled by means of p. The position and velocity gains for each mode are shown' pair- 
wise, the labelling of the abscissa referring to the mode number. Note also that a 
different ordinate scale was used for the last 50 modes to more clearly display the 
details. 

This design was obtained In less than 10 seconds computation time In the follow- 
ing manner. First, rigorous bounds for the deviation of the solution of the stochas- 
tic Rlccatl equation from the asymptotic solution are available from Ref. 3, and these 
show that the asymptotic solution Incurs negligible error (4th significant digit) for 
the 25th mode' and that the error declines rapidly for all higher modes. Consequently, 
the solution of the full-order design equation for the last 75 modes was obtained with 
virtually no computational effort. Detailed computations were required only for the 
first 25 modes. Thus, because of the special structure of the stochastic Rlccatl equa- 
tion, the results shown represent a highly accurate solution for the full 100-mode mode 


The above features Illustrate a general principle which emerges also In numerous 
additional examples: Undev a system model which incorporates limited a priori infer- 

nation, it is always possible to so arrange matters that the burden of mean-square 
optimal design computation is similarly limited. 

Also note from Fig. 10 that the design algorithm has automatically produced a 


re£isonably high gain (> 20-percent damping) , which is a virtually deterministic LQ 
design for those modes (the first three) associated with low uncertainty and an 
inherently robust rate feedback control for the poorly known modes. These regimes 
were not secured on an ad hoc basis but are limiting qualitative features of a 


global design which Is guaranteed stable and optimal . 


STOCHASTIC BEAM PROBLEM — GAIN MAGNITUDES 
(q = 0.005, p - a - 0.1) 



MODE NO. 


Figure 10 





STOCHASTIC STABILITY AND ROBUSTNESS PROPERTIES 

With regard to stability, various theoretical assurances are available for the 
maximum- entropy approach. Figure 11 sketches the basic rationale for many of the 
stability and robustness proofs. First, for a variety of controller forms, unique 
solutions have been shown to exist for all uncertainty levels. As part and parcel of 
mich demonstrations, it is necessarily shown that the mean-square optimal performance, 
J*, is bounded. But, under certain mild geometric conditions, this is Impossible 
unless, the closed-loop system is second-moment stable (and hence almost surely stable). 
In addition, one has assurance of robustness with respect to deterministic system 
stability concepts. By Lyapunov arguments and associated singular value bounds, the 
range in parameter space for which stability Is obtained can be shown to be comparable 
to the modelled uncertainty levels (as will be Illustrated presently). 

The above properties thus Imply that by allowing the modelled uncertainty levels 
to Increase without bound, one may use the maximum-entropy optimality conditions to 
discover (or rather, rediscover!) particularly robust controller, forms. As Fig. 11 
Indicates, this Is Indeed the case. Ve may summarize much of this work by stating: 
Vndep the maximum- entropy model, and in the presence of very great parameter uncer- 
tainty, the mean-square optimal control is a control which is inherently energy dissi- 
pative. 



EXAMPLE ILLUSTRATING ROBUSTNESS PROPERTIES - STABILITY STUDIES 

To Illustrate the relative stability properties of the stochastic design, con- 
sider the simple two-mass system in Fig. 12 (for complete details and further results 
see Ref. 13). With a single force actuator and displacement sensor located at mass 1 
It is desired to suppress the displacement response of mass 2 (thus the steady-state 
performance Index, J , Is as Indicated In the figure). We assume that only the fre- 
quency of the single elastic mode Is uncertain and Is subject to normally distributed 
variation with standard deviation, a, relative to its nominal value, u. Of course, 
since only one parameter Is uncertain, the range of stability In parameter space can 
be shown particularly simply for this example. 

The assumed observation noise, requires design of a full-order dynamic com- 
pensator (observer) for this system. The relevant mean-square optimality conditions 
for the maximum-entropy model are given In Ref. 13. These were solved for a range of 
values of o. Each such compensator design corresponding to a given modelled uncer- 
tainty (a given o valiie) was then Interconnected with a perturbed structural model 
with modal frequency^ u(l +6). Thus 6 Is the deviation of the frequency relative to 
the nominal value, u. Stability was then checked for the closed-loop system over a 
range of values of 6. 


TWO-MASS SYSTEM 



PERFORMANCE INDEX: 


= E[x| + u*' 

ASSUMING FREQUENCY 
IS GAUSSIAN: 

T = 

(<xS)' 


Figure 12 
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EXAMPLE PROBLEM - STABILITY BOUNDARIES 


. The results of the stability study are sunanarized in Fig. 13. Thus, for example, 
the deterministic^ LQG design (a ” 0) is unstable for relative frequency variations 
between -10 percent and -50 percent and greater than +30 percent, (to the other hand, 
progressive increase of the modelled uncertainty level dramatically widens the sta- 
bility range. In fact, as the figure makes evident, stability is guaranteed for 
d < o. Note that if 5 is randomly distributed but of bounded variation, a may always 
be chosen to guarantee stability. On the other hand, if 6 is randomly distributed 
with unbounded variation, appropriate choice of o can still reduce the probability of 
instability to any desired level. Thus, a striking Illustration of the general robust- 
ness results discussed is evident. Similar stability studies have been performed 
and analogous results achieved in other much more complicated examples (see Refs. 13-14) 


(Unstable Region Shaded) 



S RELATIVE VARIATION 


Figure 13 



MEAN-SQUARE OPTiiiAL DESIGN: IMPOSITION OF 

IMPLEMENTATION CONSTRAINTS 

The remainder of this paper considers various optimization Issues. As noted, 
we attained a maximum- entropy formulation for Implementable compensation by first 
considering the simplest controller .form, l.e., full-state feedback regulation, and 
then by progressively removing Idealized restrictions. This stage-by-stage develop- 
ment Is summarized In Fig. 14. As mentioned, the maximum -entropy model conjoined with 
mean-squat’e optimization gives rise to new forms of optimality conditions. For 
example, determination of the mean-square optimal regulator requires solution of the 
stochastic Rlccatl equation (for Its .definition and basic properties, see Refs. 3, 4, 

5 and 15). This equation reduces to the standard Rlccatl equation in the determinis- 
tic limit, but, as the foregoing beam example illustrates, new convergence properties 
are exhibited In the presence of uncertainty. 

Next, the unrealistic restriction of full-state feedback was removed and state 
Information was assumed available only through a limited number of sensor outputs. 
Assuming non-singular observation noise, the case of full-order dynamic compensation 
was considered and appropriate optimality conditions derived In Refs. 13 and 14. These 
results sho“ that in the presence of uncertainties, the separation principle no longer 
holds. Indeed, the optimality conditions determining the regulator and observer gains 
(k and F, respectively in Fig. 14) are determined by two coupled stochastic Rlccatl 
equations together with two additional Lyapunov equations. The appended Lyapunov 
equations serve to determine two auxllliary cost matrices which characterize the error 
"leaking through" from regulator to observer (and vice versa) due to the action of 
parameter uncertainty. Despite this additional complexity, efficient numerical pro- 
cedures for solution of the full optimality conditions have been developed and success- 
fully applied to a variety of example problems (of which the two-mass system considered 
above Is the simplest case). Finally, it should be noted that various desirable fea- 
tures observed for the regulator problem, l.e., stochastic simplification and the 
automatic emergence of a simple, dissipative control for the poorly known modes, are 
also obtained in this less idealized setting. 


The next and final stage In the development of suitable optimality conditions 
removes the assumption that the compensator order Is equal to that of the plant. 
Clearly, while the structure Is of infinite dimension, on-board software limitations 
force the compensator order to be rather modest. Thus It Is necessary to design a 
mean-square optimal, dynamic compensator whose order (N ,N S 2n) is preassigned by 
the designer. This Is the problem of optimal, linear, fixed-order dynamic com- 

pensation and quite apart from the Inclusion of parameter uncertainty effects; the 
results presented In Ref. 16 are new. In place of the purely computational approach 
which combines gradient calculation with parameter optimization algorithms. Ref. 16 
provided, for the first time, explicit optimality conditions for fixed-order compensa- 
tion. As in the case of full-order compensation, these conditions amount to two 
Rlccatl-llke and two Lyapunov-like equations. A qualitatively new feature, however. 

Is that these equations determine not only the optimal gains but also a projection 


(Idempotent matrix) which defines the N -dimensional observation and control subspaces 
of the compensator. This "optimal projection" essentially characterizes the geometric 
structure of a reduced-order model of the plant which is employed internally by the 
compensator. For example, the usual modal truncation approach to model reduction Is 
tantamount to the assumption that the projection is of the canonical form 
In the modal coordinate basis. However, In most Instances, the true opti- 
mal projection Is far more general than this, and, in contrast to all pre- 
vious approaches, the conditions of Ref. 16 determine both the gains and 
the optimal compensator structure with Its associated reduced-order model. 
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Various theoretical properties of the fixed-order compensation optimality condi- 
tions are elaborated in Ref. 17. In particular, solutions of these equations give 
rise to a^cojitrplj for which closed-loop stochastic stability and mean-square optimality 
are guaranteed. Moreover, solutions exist of a highly simplified asymptotic form for 
those high-order, poorly known modes constituting the incoherent range. Thus, as in 
the regulator problem, the optimization computations can be carried out for very large 
order systems. 
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OPTIMALITY CONDITIONS 

FULL STATE 
FEEDBACK 

(*A* is Uncertoin ond Treated Under 
the Max. Entropy Model) 

», 

MIN: J = C dt r {x^R.x + o^R-o) 

k “^t ' ^ 

o 

. . . «2n 

X •• Ax + Bo + Wj . X c R 

i “ — kx . o f R* 

/ 

• SPECIAL FEATURES 
STOCHASTIC RICCATI EQ'S 

• NEW CONVERGENCE 
PROPERTIES WITH 
INCREASING MODEL ORDER 

FULL-ORDER 

COMPENSATION 

j'vME AS ABOVE, EXCEPT: 

u = — kq . o « R* 

. , .2n 

q = oq + Fy . qrR 

p 

y = Cx+W 2 ; yCR 

COUPLED STOCHASTIC 
RICCATI EQ'S 

• NO SEPARATION PRINCIPLE 

FIXED-ORDER 

COMPENSATION 

SAME AS ABOVE, EXCEPT 
N 

q f R ** , N^ ^ 2n 

COMPENSATOR STRUCTURE 
DETERMINED Bf OPTIMAL 
PROJECTION 


Figure 14 










»• 


•t 


I ' 



91 


f ' 


. fi > 





GENERALITY OF THE OPTIMAL PROJECTION APPROACH TO FIXED-ORDER COMPENSATION 

A further property of the formulation of Ref. 16 Is Its considerable generality, 
as Indicated in Fig. 15. Although the original derivation assumed nonsingular obser- 
vation noise and thus permitted no direct sensor feedback, the conditions of mean- 
square optimal output feedback can be recovered from the equations of Ref. 16 as a 
special limiting case. Also if N is set equal to the plant dimension, the optimality 
conditions immediately reduce to ?hose of Ref. 13. If, in addition, no parameter 
uncertainty exists, the full-order compensation optimality conditions directly yield 
the uncoupled Riccati equations of LQG theory. Finally, a variety of LQG-based tech- 
niques with their assorted model reduction schemes can be seen as special approxima- 
tions to the present optimality conditions (see Ref. 17 fbr further details regarding 
the relations between the optimal projection and mdern nx>dal control approaches). 
Thus, the optimal projection formulation has the conceptual scope needed to assess 
the relative merits of the various more specialized design schemes that have been 
advanced in the recent past. 





MAXIMUM-ENTROPY MODELLING CONJOINED WITH OPTIMAL PROJECTION; 

TOWARD MEaiANIZED DESIGN SYNTHESIS 

Although numerical techniques for solution of the fixed-order compensator opti- 
mality conditions are based primarily upon software previously developed for full- 
order compensation. It remains to consolidate this software and explore fixed-order 
design for representative example problems. The theoretical properties mentioned 
certainly provide strong Impetus for completing this task, and our Immediate object 
Is the design synthesis software package indicated In Fig. 16. Given a high-order 
plant model (specified nominal values and uncertainty relaxation times), together with 
quadratic control objectives and a specified compensator order, this software package 
would handle the optimization computations via a two-level Iterative method employing 
efficient relaxation techniques for the Incoherent regime. Because of the guaranteed 
system properties listed at the bottom of Fig. 16, the need for direct human interven- 
tion In the control design process would be largely eliminated. The designer's 
burden is thereby greatly reduced, and all that Is required on his part Is a familiarity 
with the meaning and Impact of uncertainty relaxation times and some design experience 
with quadratic optimization. 
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OBJECTIVES 
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Figure 16 
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UNIFIED STRUCTURE/ CONTBOl. DESIGN SYNTHESIS 

The contemplated design optimization package would be of particular usefulness 
In performing trade-off studies needed for the overall structure/control design syn- 
thesis (see Fig. 17). In the very earliest design stages, the theory underlying the 
control design package can be used to delimit the various qualitative regimes and 
bound the best achievable nean-square performance for several candidate structural 
design concepts. This would allow computationally simple performance trade-off 
studies- including the effects of a pnori modelling uncertainties, actuator/sensor 
placement, etc. Such studies would permit the early elimination of those structural 
concepts which are not feasible in view of the performance requirements. In general, 
the control design approach outlined here would help to guide the structural design 
according to the Inherent capabilities and limitations- of fixed-order compensation. 

In the later stages of design (with a reasonably well-defined structural concept), 
the control synthesis soft. rare could be used to achieve a dynamic controller of lowest 
order which meets performance requirements. If the resulting compensator order also 
satisfies on-board computational constraints, a final linear system design is achieved 
which produces acceptable performance despite parameter uncertainties taken at their 
a priovi levels. In this case, the need for parameter Identification twuld be obvi- 
ated. Generally, however, a priori, uncertainties may be so large that the robustness 
level produced by the control design synthesis package Is too high to permit good 
mean-square performance. Here, the design synthesis package of Fig. 16 can be 
employed Iteratively to establish the gradients of performance with respect to the 
various relaxation times. This would permit one to enumerate those (relatively few) 
"crucial" structural parameters that must be identified In order to secure the 
required performance. Considerable simplification thereby results for Identification 
algorithms, since only the crucial parameters need be addressed. 
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Kinetic Energy: T(t) ■ y m(P)u^(P,t)u(P,t)dD(P) 

u(P,t) ■ displacement vector 
m(P) ■• mass density 

D ” domain of extension of the structure 
P • nominal point 

Potential Energy: V(t) ■ ^ [u, u] 

I , ] “ energy inner product (Ref. 1) 

Virtual Work: 6 W - f'’’(P,t) 6 u(P,t) dD(P) 

f(P»t) “ distributed force (actually explicit function of state rather than of 
spatial position and time) 
du(P,t) " virtual displacement 


Extended Hamilton's Principle 


■i:: 


(«L + 6 H)dt - 0, 

6 u(P.t) « 0, t “ tj^, t^; P 6 D 


L ■ T - V * Lagranglan 

Partial Differential Equation of Motion for the Structure (Ref. l)i 

iu + Mu - f , P 6 D 

L » Differential operator matrix with entries of order 2p 
M ■ Mass matrix 

Boundary Conditions: B^u “ Q , P 6 S; 1=1, 2, ...,p 


B 


^ - Differential operator matrices with entries of maximum order 2p-l 


S " Boundary of D 


STRUCTURE DISCRETIZATION BY THE FINITE ELEMENT METHOD 

Linear Transformation: u(P,t) = L(P) 5 (t) 

L(P) “ matrix of Interpolation functions (Ref. 1) 

3(0 ” vector of nodal coordinates 

Discretized Kinetic and Potential Energy: T ■ ^ q^Mq » ^ ^ q^Kq 

M - /jj L^^(P) ML(P) dD(P) - mass matrix, M - m’^ , M>0 
K “ (L(P) , L(P)1 “ stiffness matrix, K ” k!^, K^O 

T 

Discretized Virtual Work: 6 W « Q 6 q 
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9 ** /d ^ (P)f (P,t)dD(P) “ nodal force vector 
Discretized Equations of Motion: Hq + Kq ~ Q . 

Equations have appearance of an n-degree-of -freedom discrete system. 

NATURE OF THE DISCRETIZED SYSTEM 

A discretized model Is a truncated model meant to represent a distributed structure. 
Eigenvalue Problem of Order n: , r*l,2,...,n 

K' 


" n X n stiffness and mass matrices 
X^“^ - computed eigenvalues, * * •— 

n-dlmenslonal eigenvectors 


3 ^ 


j^(n+l) (n+1) (n+1) 

r -r 


Eigenvalue Problem of Order n + 1: 

r — 1,2,... ,n+l 
.(n+l)^-(n+l)^ 

*1 -2 -••-Vl 

Incluslc.; Principle (Ref. 1): 

X — 1 — 2 — 2 — — n — n — n+1 

Convergence Property: 11m X^°^ - X 

n+« ^ 


r* 1,2,. a a ,n 


Paradox: By Increasing n, more computed eigenvalues tend to be accurate, but new 

ones are added at the high end of the spectrum, and the latter ones tend 
to be wildly In error. 

Conclusion: No discretized model can yield a complete and accurate representation of 

a distributed structure. 

MODAL EQUATIONS OF MOTION FOR CONTROL 

Consider an n-degree-of-freedom discretized model, but drop superscript (n) , 

X^ “ computed eigenvalues, q^ - computed eigenvectors 
Computed Eigenfunctions for Distributed Structure: <|>^(P) ” L(P)q^, r**l,2,...,n 

Structure is self-adjoint — >■ modes are orthogonal (Ref. 1) and can be normalized 

9l«Sr**rs* «„• f.8 - 1.2,...,n 

Consequence: /jj ♦g(P)M*^(P) dD - dD - r,b - 1,2 n 


1 
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H . 

Expansion Theoren (Ref. 1): u(P,t) “ X ♦ (P)u (t) 

.r-1 ■■ ■" . 

o^(t) ■■ modal coordinates 

N “ number of modes retained In the model, N«n 

Hodal Equations: u(t) -f Au(t) “ f(t) 

T 

5 (t) ” [uj^(t) UjCt) ... Ujj(t)] “ modal displacement vector 

f(t) ” [fj^(t) ” nodal control vector 

2 2 2 

A - dlag (Xj^ Xj... Xjj] » dlag ••• “jjl ” matrix of eigenvalues 

* natural frequencies 

Modal Controls: f^(t) - *^(P)f(P,t) dD, r-l,2,...,N 


M 


COUPLED MODAL CONTROLS 
Feedback Control Vector: f ” f(u, u) 

Distribut'd Controls In Terms of Discrete Actuators: f(p,t) ■ J F. (t) 6(P-P.) 

J-1 J 

F. (t) “ discrete-point controls, i<»l,2,...,M 

-J t. o 

M T “ T 

Modal Controls: f (t) - I L *,.<P)F. (t)6(P-P.)dD - I *^(P.)F. (t) , r-1,2 , . . . 

u j -r j 

Control Spillover (Ref. 2) Into Unmodeled Modes: 0, r“N+l ,N+2 , . . . 

Modal Feedback Control: f •• B F 


B ■■ Modal Participation Matrix, B 


Modal Equations: u(t) + Au(t) * B F(t) 


. - [$^(Pj^) $^(P2)... $^(P„)] 


M modeled nodes 


controlled modes 


' residual modes 


, N-N^ + Nj^ 


Hodal Equation In Partitioned Form: 


Feedback Control Vector: F “ Gj^ 

C^, “ control gain matrices 


.h. 

+ 

o 

1 

1 — 

1 

Si I 
( ! 

m 

1 

Ol 
M ! 

1 

a 


r i** J 

L'^rJ 

1 

w| 
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f 
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Clon>;d-Loop Modal Equations 




Closed-loop poles of the controlled modes are determined by -B G, and A -B G . 

C1CC2 

Closed -loop poles of the residual modes are determined by 

Conclusion: Closed-loop poleti of thie residual modes are not affected by feed- 

back controls; therefore, control spillover Into the residual 
' ‘ \ modes cannot destabllze these modes. 


CONTROL IMPLEMENTATION USING OBSERVERS 
State Equations of Motion: 




Feedback Control Vector: F 

•T T,T 






yc ” fi?c i?c^ 


modal state vector 


’Oj-A^ ■ 

L 1 1 0 j 






®c ®R 

• ®c “ ' K “ 

^0 0 


Displacement Measurements Vector: y(t) - Cu(t) * C^u^Ct) + Cj^Uj^(t) 

Velocities Measurements Vector: y(t) - Cu(t) - Cj,Ug(t) + C^u^(t) 


*(P2> 

♦ (V 


♦ (P^) 


modal matrix evaluated at P 
K ** number of sensors 


Output Vector: y(t) - C*v^ + C*v 

The term C^Vj^ gives rise to observation spillover (Ref. 2) 

[ Yc T r I ° ^c ®c 1 

yr J L® JlyrJ [k I 


_^c_ 


(y-?) , 


y • ^cYc + 


R 
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Vg, Yj^ - estimates of y^, Vj^ 

K^, ^ * observer gain natrlces^ chosen so thbt approaches expor>sntlaliy 
y “ observer output 


Feedback Control Vector; F 


Vc 


- v« - v„ 


Error Vectors; * Yc “ Yg » Sr ir 
S tate and Error Equations Combined; 


Yc 


A 

jc. 


Ac + ^c^^c 1 

0 ; 

®c°c 1 0 

Vc i 


i « 

0 i 

1 

0 j 

Ac “ *^C^C 1 ~^cS 

1 

0 1 

0 j 

-V‘ I \-hPk 


-I5. 

-:L 

.!s_ 

?R 


Conclusion: Observation spillover exists but cannot destabilize the system If the 

observator gains are chosen properly. 


fi-: 


9 * 

li 


CO'JVROLLABILITY AND OBSERVABILITY 

! A^B; I I B^l 

For system to be completely controllable, C must be of rank 2Np(Ref. 3). 


Controllability Matrix: C ■ [B^ | A^B^ J A^.B^ j , ..j, 


Controllability Matrix for Discretized Model: 

i I 

I 


“c i° iVc 


I 

0 ; B 


C i ° 


-"c»c I 0 


(_1)NC-1 ANc -1 


For system to be completely controllable, all rows of B^ must be nonzero 
Observability Matrix: 0 


[C; i aX c'l I (aJ)2 c'^ i (aA) c'Al 


.T.^”c-l _,T, 


,T { T 

J '^c ^ c j '"c' " C I '"C' “ C‘ 

For system to be completely observable, 0 must ha .’e rank 2N^(Ref. 3), 


Observability Matrix for Discretized Model: 

0 



1 ^ 

i -*c 

1 a 2 c’' 1 
1 C ‘'c 1 __ 

0 

0 

j -vl 

1 ” 

a" 1 0 1 

a”c-ic! 


C C 1 j 

c c 


For system to be completely observable, all columns of C^ must be nonzero. 
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MODAL FILTERS (IN SPACE) 

Second Part of Expansion Theorem: u (t) ■ /^♦^(P)Wu(P,t) dD 

' ~ r-1,2,... 

u(P,t), u(P,t) ” distributed measurements 

Integrals filter out contributions from modes r, hence, the term modal filters 
(Ref. A). 

Dli.’crete Measurements: Use K sensors to produce discrete measurements. Then, use 

the finite element method to generate estimates 
u(P,t), ,6(P,t) of u(P,t), u(P,t). 

Modal Filters for Discrete Measurements: Uy(t) •• /p^^(P)Mu(P,t) dD 

A A T * r“l,2,...,N„ 

Ur(t) - /jj$J(P)Mu(P,t) dD 

Finite Element Approximation: Uj^(P,t) ■ L^(P)y^(t), Uj^(P,t) ■ b'*'(P)yj^(t) , P 6 D^^ 

1 - 1,2,...,K-1 
“ nodal measurements vectors 


Modal Filters: u (t) - I L <|>^(P)Mu, (P,t) dD - ^ I y.(t), u (t) =* I I y (t) 

I^^ - /p ^^(P)ML^(P) dD (computed off-line) 

INDEPENDENT MODAL SPACE CONTROL (IMSC) METHOD 

(1) Control of Distributed Structures by Distributed Actuators 

2 

Infinite Set of Modal Equations: u^. - r - 1,2,... 

Special Case of Modal Feedback Controls: f^ - fj.(“^» ^ 

Note: The rth modal control depends only on the rth modal coordinate and 

velocity. 

Consequence: System becomes c. .pletely decoupled ; therefore, the control for 

each mode can be designed independently of any other mode. 

Synthesis of Actual Distributed Controls from ffod. 1 Controls (Ref. 4): 

f(P.t) - I M(P)d (P)f (t) 

■ r-1 . ■ 

Note: If only Nq modes are to be controlled, then set f^ - 0 (r > N^;) , which 

implies no control spillover Into the uncontrolled modes. 






I 


■*A 


: ''.'r V 


“ control [jaln parameter 


Linear Controls: ®r2**r* *■ “ 

0^2 “ modal gains ‘ . 

Hbnllnear, On-Off Controls (Ref. 5): “ ^“^r* ^r ^ **r* l“jl ^ ^ 

r ~ 1,2,. . . 

2d^ ■ magnitude of dead-band region 
(11) Control of Distributed Structure by Discrete Actuators 
Equations for Control of Modes: “ ®C- 

Synthesls of Actual Controls from Modal Controls: F “ ®C~C * ^ ^ **C 

t - ■ «C 

In the latter case, the number of actuators must be the same as the number 
of controlled modes. 

•• —1 

Equations for Uncontrolled Modes: Uy + - ByF “ ®jj®q Cq 

Control spillover Into the uncontrolled modes does exist. 

No Instability Is possible. 

(Ill) Control of Discretized Structures by Distributed Actuators 

«C 

Sase as (±) , except, that f(P^^;' “ J M(P)$ (P)f (c) 

r-1 ^ 

No control spillover Intc the uncontrolled modes exists. 

(iv) Control of Discretized Structures by Discrete Actuators 
Same as (11) . 


COMPUTATIONAL ALGORITHMS FOR CONTROL 

(1) Pole Allocation 

Select Closed-Loop Poles pj^ (fc-l,2,...,2Ny) 

State Equation for Single Input: v_ ■ A_v„ + b' F(t) 

C*’C "'C 

Single Input Feedback: F(t) “ 

gy ■ gain vector 
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2Hq 

klA-Vsj 

^ ' J-l I ^"c 


Yj , . . . ,2Ng) ” left eigenvectors of 

General Multi-Input Control Is Not Feaeible. 

Dyadic Control (Inputs Are Proportional to One Another) : F(t) » G v (t) 

^ c**c 

■ gain matrix 


Gc-h Z 


!C k-1 


ZNc - 

* (Pt - *^)v. 


1«1 T 

J c J^l k J 


Closed-Loop Poles for IMSC; 


Pr - t*r - *®r* 


Modal Controls: f^(t) - - a2 . + 2a^u^. r-1,2 


(11) Linear Optimal Control 


P®i^^ormaPce Index for Coupled Controls: J “ f ^(v^Qv + F^RF) dt 

m *0 "G 'G ' - 

Q. R ■ Weighting matrices 

>. Optimal Control Vector: F(t) - ~R*'^''*'K(t)v^(t) 

^ K(t) - 2Ng X 2N^ matrix satisfying the Rlccatl equation 

K ” -KA^ - A^K - Q - KB'R"^B’\ 

Steady-Stace Case: K - 0 

Nonlinear algebraic matrix Rlccatl equation can be transformed to eicen- 
■? solution of a real general matrix. 

• 3 

^^8o*^l-thm requires 600 multiplications for convergence. 

Transient Case: K 0 ' 

f Nonltoear matrix Rlccatl equation can be transformed into a 4Np x 4Np linear 

^ matrix equation, which must be integrated on-line. Process requires ample 

computer capacity, %dilch may make real-time implementation Impossible. 





Perfornaiice Index for THX: J “ 


f^r 

r»l 


« independent perfortnance Indices 

J, - Yr ^ V^r^***"* 


V 

»C 


<J^ — 2 X 2 -diagonal matrix, - scalar 

Solve N„ Katrlx Rlccatl Equations of Order 2. 

c 

Steady-State Case Requires N^/2 Operations, Fewer by a Factor of 1200 Than 
Coupled Controls . 

Transient Case Implementation In Real Time Is Possible. 


COMPARISON OF COUPLED CONTROLS AND IMSC 

Advantages of Coupled Controls: Fewer actuators, provided controllability Is 

ensured. 

Advantages of IMSC: 

a) Larger choice of control techniques. Including nonlinear control 

b) Lower computational effort (Table 1) 

c) Lower computer storage requirement 

d) Lower control energy (Ref. 6) 

e) Locations of actuators not Important (Ref. 7) 

f) Provable robustness (Ref. 8) 

Some of these advantages are Illustrated In Figures 1-5. 


Number of 
Actuators 

Time 

(b) 

.3 

0 

4 

176.5 

M 

u 

C 

6 

127.1 

0 

0 

8 

104.5 

•t> 

10 

100.0 

tu 

9 

0 

0 

12 

87,7 

IMSC 

12 

2.3 


Table 1.- Compatatlonal effort required to solve the Rlccatl equation. 
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1 actuator 
5 actuators (IMSC) 


Figure 1.- Displacement at the end of the bar - pole allocation 


1 actuator 
5 actuators (IMSC) 


KiTtrii 


Figure 2.- Energy required for the controlled modes - pole allocation 


6 actuators 


4 actuators 


6 actuators 


in actuators 
12 actuators (IMSC) 


« o.uu iu.wj isj.ga zo.oo Z5.00 30.00 

T - ■ 

Figure 3.- Performance Index - coupled and Independent controls. 
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ANALYSIS AND TESTING 
: SPACE STRUCTURES 


C. V. Stable 

Space Systems Division, General Electric Company 
Valley Forge, Pennsylvania 



GROUND TEST CONSIDERATIONS 



Figure 1 outlines key considerations for performing ground tests of LSS. The large 
size combined with the loading due to gravity will make testing of the complete struc- 
ture difficult. Gravitational stiffening, suspension effects, virtual air mass, pre- 
loads, and air damping will alter the dynamic characteristics. Low resonant fre- 
quencies and high modal densities within the frequency range of Interest combine with 
small motions and accelerations to make testing difficult. Mechanism complexities 
and nonllnearltles associated with space-erected /assembled structures cause struc- 
tural complexity regardless of other considerations. From these considerations it ap- 
pears that ground testing of a complete LSS will be difficult if not Impossible. 





LARGE SIZE 

GRAVITY - STIFFENING, LOADING 
AIR MASS, DAMPING 
LOW RESONANT FREQUENCIES 
HIGH MODAL DENSITY 
SMALL MOTIONS/ACCELERATIONS 
SUSPENSION EFFECTS 
MECHANISM COMPLEXITIES 
NON LINEARITIES 


GROUND TESTING WILL BE VERY 
DIFFICULT IF NOT IMPOSSIBLE 


Figure 1 
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Figure 2 outlines some approaches to ground testing. Scale model testing can ever- 
come the size and gravity effects but does not assure detailed representation. Element 
tests provide a means of supplementing model tests and will be particularly effective 
for Items such as Joints. Substructure testing appears to be the most promising, en- 
abling dynamic behavior of actual structure to be measured while overcoming some of 
the concerns (e.g. , size, gravity, modal density). Intentional linearization of the 
structure will be helpful In minimizing test problems associated with small orbital 
motions. A separate analysis of the test configuration will probably be required in 
That some ground teat effects will undoubtedly be present. Most likely, tests will be 
'.irected toward verification of analysis using either models or partial structures. 


• PERFORM SCALE MODEL TESTS 

• PERFORM ELEMENT TESTS 

• TEST SUBSTRUCTURES 

• LINEARIZE STRUCTURES 

• PERFORM SEPARATE ANALYSIS OF GROUND TEST ARRANGEMENTS 


TESTS WILL BE DIRECTED TOWARD 
ANALYTICAL VERIFICATION USING 
MODEL OR PARTIAL STRUCTURES 


Figure 2 
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ROLL-DF SOLAR ARRAY MODAL TEST 



The modal test of a roll-up solar array performed at GE, Figure 3> indicates test ap- 
proaches which have been used. The test was performed in the large vacuum chamber 
(32 feet diameter by 54 feet high). The array was hung downward to eliminate suspen- 
sions and preclude excessive gravity loading. Sinusoidal base motion was used to 
excite the array. Special low-frequency test equipment was used including an oscil- 
lator/analyzer that would provide a "hj’perbollc" sweep (log rate Increasing with fre- 
quency to provide constant resonant sweep distortion) . Air caused a 60-percent re- 
duction in the fundamental frequency and an order of magnitude change in damping. 
Gravity, Included in the analysis, shifted the fundamental resonant frequency by 
a factor of 2. Low-frequency sinusoidal testing was found to be tedious and very 
time consuming. 



• DC COUPLED MICH-DISPLACEMENT VACUUM RA 
VIBRATION EXCITER 

• LO«H=REQUENCY SWEEP OSCILLATOR 

•oet Hz 

LINEAR, LOGARITHMIC OR HYPERBOLIC SWI 

• NONCONTACTINC VACUUM RATED SENSORS 

• CO/QUAO RESPONSE MEASUREMENT 

• BASE MOTION EXCITATION 

• NO SUSPENSION EFFECTS 

• ANALYSIS INCLUDED GRAVITY 

• AERODYNAMIC EFFECTS SHIFTED THE FUNDAM' 
RESONANCE «0 PERCENT 

• GRAVITY SHIFTED THE FUNDAMENTAL RESON 
BY A FACTOR OF 2 


TEST ARRANGEMENT 


Figure 3 
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DESIGN DEVELOPMENT AND VERIFICATION 


The development and verification process for LSS will differ from that currently being 
used for spacecraft as indicated in Figure 4. Present practices rely primarily on 
ground tests to verify the design and verify workmanship for the flight unit using a 
prototype (qual) structure. Testing of the complete LSS will not be practical. Ground 
tests will be used primarily to verify analysis, and the use of ai prototype structure 
will be eliminated to reduce cost. A key feature will be that the ground test verifi- 
cation role will be replaced by analysis. Consequently, major emphasis must be placed 
on the development of accurate analytical methods. 



Figure 4 
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ANALYSIS CONSIDERATIONS 



The analytical considerations can be grouped in the areas of structural Interactions 
and complexity as Indicated In Figure 5. The structural Interactions can involve 
motion-dependent forces caused by orbital rates » the control system, and solar radia- 
tion. All of these interactions have been experienced by spacecraft. Spacecraft with 
highly flexible gravity gradient rods exhibited destructive orbital Interaction, and 
the performance of other spacecraft was degraded by thermal "flutter." Structural 
analysis complexities Include FEM modeling of distributed mass systems requiring large 
DOF or condensation, modeling adequacy for mechanism Joints, inaccuracies arising from 
normal mode analysis of structures having nonproportional damping, and analysis of non- 
linear structures such as those using tensioning members. Comprehensive LSS analysis 
will require complex models and use modified equations of motion that Include Inter- 
actions. 





• STRUCTURAL INTERACTIONS 

- ORBITAL 

- CONTROLS 

O 

- THERMAL 

• STRUCTURAL COMPLEXITY 

- DISTRIBUTED MASS 

- MECHANISM MODELING 

- DAMPING 

- NONLINEARITIES 
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. ANALYTICAL APPROACHES 



Analytical approaches which could be used are outlined in Figure 6. To evaluate in- 
teractions, a combination of computer codes that treat each structural interaction In 
detail and a comprehensive code that treats all Interactions simultaneously appears, 
to provide a tractable approach. Specified codes would be used to support the struc- 
tural design and would be more detailed, more readily Implemented, and less expensive 
to run. The comprehensive analysis would build on the specialized design support 
analyses and be used to verify the design. Attractive approaches to modeling of the 
complex LSS configurations Include substructuring, modal synthesis, and experimentally 
based models. Modal models are attractive in that they provide frequency truncation. 
An Integrated modeling approach that considers experimental verification should be 
used to maximize confidence in the analytical predictions. 



j; • DEVELOP STRUCTURAL INTERACTION CODES 

^ - ORBITAL 

- THERMAL 

- CONTROL 

- COMBINED 

• DEVELOP COMPLEX MODELS OF STRUCTURE 

- SUBSTRUCTURE ANALYSIS 

- MODAL SYNTHESIS 

- EXPERIMENTALLY BASED MODELS 






IWTOEM SPACECRAFT ANALYSIS ^ 

The analysis and test approach being used for a tandem-launched spacecraft is shown in 
Figure 7 and illustrates a substructuring approach. The spacecraft are launched tiro at 
a time with the one Inverted relative to the other. The analytical model is assembled 
from six substructures: two identical main bodies, four identical solar array stacks. 
The modal test will provide solar array modes and main body modes that can be used to 
synthesize the tandem spacecraft pair. One solar array stack was tested and two main 
body tests will be performed. Residual mass and flexibility from the FEM will be used 
with the test erodes to assemble the spacecraft pair. Tests of a single solar array 
and main body simplify the modal tests and reduce the amount of test hardware required. 


LAUNCH CONFIGURATION 
SOLAR ARRAY 



irKZXZJTI 


L 

MAIN 

BODY 


s: 


Is 


MAIN 

BODY 



ZTA 

1: 


UPPER 

S/C 


LOWER 

S/C 


MODAL SYNTHESIS 
ANALYSIS 


ncxin 


MODAL TEST 
ARRANGEMENTS 


15ZSZT 

7 ^ 7 ^ 


ll 


UUZILl 

7 ^ 77777777 ^ 


STRUCTURE SYNTHESIZED 
FROM 6 SUBSTRUCTURES 



STRUCTURE TESTED 
USING TWO SUBSTRUCTURES 


POST TEST ANALYSIS USES 
MEASURED MODES WITH 
RESIDUAL MASS AND 
STIFFNESS FROM FEM 


Figure 7 








Some research areas are outlined In Figure 8. The first area Is the development of 
analytical codes that will enable accurate prediction of orbital dynamics behavior In- 
cluding structural Interactions with orbital, thermal, and control excitations. As In- 
dicated previously, specialized codes and a comprehensive code are suggested to sup- 
port spacecraft design and verify orbital performance. Ground and orbital test sub- 
stantiation of the codes Is an essential part of code development. Much of this 
activity Is currently planned or Is in progress. 
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Because of the Increased reliance on analysis. It Is essential that testing and methods 
of using experimental data to enhance analysis accuracy be developed. Also, because 
conq>lete structures will probably not be ground tested, the development of methods for 
^ using data generated with partial structures is essential. Specific items needing 
4 research are (1) methods of using modal test data to Improve FEM models, (2) methods 
of correcting test data to eliminate ground test effects, (3) limits for ground tests 
which assure adequate data for correlation with analysis, (4) methods of performing 
modal tests (ground or orbital) vyhen modal density Is high, and (5) methods of analysis 
and testing using substructures. Damping should be examined as a means of reducing 
structural dynamic sensitivity. 


• ANALYTICAL CODE DEVELOPMENT 

- COMPREHENSIVE (STRUCTURE. ORBITAL, CONTROL. THERMAL). 

- SPECIALIZED (STRUCTURE-ORBITAL, STRUCTURE-CONTROL, "ETC.) 

- TEST VERIFIED (GROUND, ORBITAL) 

• ANALYSIS/TEST CORRELATION METHODS 

• SUBSTRUCTURE ANALYSIS/TEST METHODS 

• TEST METHODS 
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THERMAL ANALYSIS CONSIDERATIONS FOR 
LARGE SPACE STRUCTURES 


Howard M. Adelman and Charles P. Shore 
■ NASA Langley Research Center 
Hampton, Virginia 


THERMAL-STRUCTURAL ANALYSIS CONSIDERATIONS FOR VARIOUS MISSIONS 

Some of the aerospace vehicles which have been the subject of NASA research are 
shown In figure 1. Each vehicle Is associated with a different set of thermal 
loading circumstances. These circumstances determine the need for thermal -structural 
analysis tools. The supersonic transport concept presented minimal analytical 
difficulties* and it was found that existing analysis tools were adequate. The 
hypersonic transport design included an actively cooled structure to dissipate high 
heat loads of long duration. This configuration required analyses to account for 
strong coupling amcng thermal, structural, and fluid behavior. Studies of this 
vehicle led to the development of an integrated thermal -structural analysis approach 
and integrated finite elements in which the temperatures and thermal input to the 
structural analysis are computed together in a consistent iianner (ref. 1). Thermal- 
structural analysis cf the Space Shuttle orbiter during entry was probably the most _ 
severe challenge to thermal -structural analysis techniques to date. Shuttle analysis 
requirements led to research to implement faster transient techniques, improved 
modeling methods, and use of advanced computer hardware (ref. 2). For large space 
structures, the impact on thermal structural analysis methods appears to be a need 
to efficiently handle radiation heat transfer effects Including view factors, solar 
flux, and shadowing. Additionally, prediction and control of the thermal deformations 
of sensitive components such as antennas will influence both analysis and optimization 
methods. 


SUPERSONIC TRANSPORT 



ANALYSIS TOOLS ADEQUATE 


HYPERSONIC TRANSPORT 


COOLED STRUCTURE - THERMAL /FLUID/STRUCTURE 
HIGHLY COUPLED LEO TO INTEGRATED 
FINITE ELEMENT ANALYSES 


SHUTTLE 



^ analysis STRAINED METHODS 
WTERPOLATION USED LEO TO RESEARCH M FASTER 
ALGORITHMS. USE OF ADVANCED COMPUTERS. ETC. 



. > ; RADIATION DOMMATEO SHAPE CONTROL NEEDED 


Figure 1 






Much of the recent Langley research In thermal -structural analysis methods was 
motivated by problems encountered In the design of the Space Shuttle orbiter (fig. 

2) for which, it was necessary to have knowledge of the histories of temperature and 
thermal stress due to the time-dependent reentry aerodynamic heating. The straight- 
forward approach would be to use a finite element and/or lumped-parameter model of 
large sections or complete components of the vehicle structure. However, due to the 
shear size of the required model and the high cost and long computer run times. It 
was necessary to model small portions of the structure. Indicated In black, and 
interpolate results to obtain temperatures In the unmodeled (white) portions of the 
orblter. About 90 percent of the temperatures used In stress calculations were 
obtained by Interpolation. There was some uncertainty and uneasiness associated 
with the use of Interpolation. Examination of analysis requirements Indicated that 
certain complicating factors listed on the figure were primarily responsible for the 
excessive resources requirea for the analysis. Research described In reference 2 
was motivated by the need to reduce the resources associated with the cited factors. 
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ASPECTS OF THERMAL ANALYSIS OF LARGE SPACE STRUCTURES 

A$ a guide to subsequent discussions of research opportunities In thermal- 
structural analysis of large space structures, figure 3 depicts various aspects of 
the analysis. Key areas of analytical deficiencies Include efficient view factor 
and flux calculations and proper accounting for Interactions. Two Important 
Interactions cited are the need to simultaneously monitor and control the 
temperature distribution and shape of flexible orbiting structures such as antennas 
and the need to account for the coupling between thermal deformations and radiation 
view factors In flexible structures. 
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MTIVATING QUESTIONS AND ISSUES 

Previous experience In thermal -^structural analysis activities for aircraft and 
space transportation vehicle structures suggests that In thermal -structural analyris 
of large space structures (Including space stations), a number of critical questions 
and Issues tend to guide research. As Indicated In figure 4, computer time and cost 
are Important factors. Consequently, Improved solution methods (algorithms) and use 
of simplified analysis and modeling techniques, such as Isothermal finite elements, 
are of Interest. The roles of advanced computer hardware, such as Interactive 
minicomputers and super mainframes, need to be assessed by comparative analyses. 

Two Issues are cited: (1) the need exists to demonstrate new promising analysis tool 

on real problems of sufficient size and complexity to ensure that good performance is 
not restricted to small academic problems; and (2) a certain dichotomy exists, 
namely, much ot the research Is carried out in the Integrated thermal -structural 
finite element context while much practical thermal analysis uses lumped-parameter 
methods. Significantly more dialogue between the two communities Is needed to 
ensure that needed and usable research Is being carried out. 


QUESTIONS 

HOW CAN WE OVERCOME EXCESSIVE REQUIREMENTS OF COMPUTER 
TIME AND COST FOR THERMAL ANALYSIS OF LSS? 

HOW DETAILED MUST ANALYSES AND MODELING BE FOR ACCEPTABLE 
RESULTS? 

WHAT IS THE ROLE OF ADVANCED COMPUTER HARDWARE? 

ISSUES 

NEED TO DEMONSTRATE TECHNIQUES ON REAL PROBLEMS IN ORDER 
TO ASSURE CREDIBILITY. 

NEED FINITE ELEMENT - LUMPED-PARAMETER DIALOGUE 


Figure 4 



;; MOTIVATION FOR USING FINITE ELEMENTS 

Recent thermal -structural research activities at Langley have been based on 
finite element analyses (fig. 5). The primary consideration for using finite 
element methods has been the need for Integrated thermal -structural analysis and 
optimization of heated structures. Use of finite elements permits a high degree of 
compatibility between thermal and structural models and avoids the cumbersome data 
transfer often required between lumped- parameter grids for temperatures and finite 
element grids for thermal stress analyses. Further, the graphics available in 
finite element codes permit rapid model checking and verification. The SPAR finite 
element computer program Is used extensively at Langley for both thermal and 
structural analyses. It Is a production level program with a modular configuration 
which gives a high degree of flexibility in terms of processor execution sequence. 

The modularity and the flexibility also enhance the program for interactive computing. 
Because of these features, SPAR is utilized as a test bed for implementation of new 
methods at Langley. 


• THERMAL-STRUCTURAL INTEGRATED ANALYSES 

• COMPATIBILITY 

• EASE OF DATA TRANSFER 

• MCDELVERIFICATION AND CHECKING 


SPAR PROGRArA USED AT LANGLEY 

• INTERACTIVE MODE 

• FLEXIBLE AND MODULAR 

• TEST-BED FOR NEW METHODS 

• PRODUCTION LEVEL PROGRAM 
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The Langley approach to improved thermal -structural analysis methodology 
consists of several parts (fig. 6). To begin, needs are identified from various 
sources including industry and acadeMic contacts and from working aircraft and 
spacecraft, problems. Once a need has been identified, ideas are conceived to 
satisfy the need, and the methodology is developed to implement these ideas. The 
methodology Is then evaluated on small well-defined problems in study computer 
codes. After the methodology reaches sufficient maturity, it is installed in SPAR 
and exercised on sufficiently complicated problems to demonstrate the benefits 
derived from the new methodology. An example of this approach is the recent 
addition to SPAR of the GEARIB set of solution algorithms as described in reference 
2. When the new methodology is independent of the analysis method (i.e., finite 
element or lumped parameter) it should also be useful in lumped -parameter programs, 
such as MITAS or SINDA. 


OVERALL APPROACH 
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ELEMENTS OF RESEARCH PROGRAM FOR THERMAL ANALYSIS AND OPTIMIZATION 

The list of research tasks shown In figure 7, apart from the division into ; 
analysis and optimization, falls into three categories: (1) items which are generic 

and ongoing (solution methods for transient temperatures, integrated analysis 
methods, evaluation of advanced computer hardware, sensitivity analysis, and 
thermal -structural optimization); (2) items which are ongoing but will receive 
Increased attention because of their importance to large space structures 
(approximate analysis techniques, improved radiation analysis); and (3) new areas of 
work (finite element modeling and analysis techniques for heat pipes). 


ANALYSIS 

• SOLUTION METHODS FOR TRANSIENT TEMPERATURES 

• INTEGRATED ANALYSIS METHODS 

• APPROXIMATE ANALYSIS TECHNIQUES 

• IMPROVED RADIATION ANALYSIS 

c 

• EVALUATION OF ADVANCED COMPUFER HARDWARE 

• MODELING TECHNIQUES FOR HEAT PIPES 

OPTIMIZATION 

• SENSITIVITY ANALYSIS 

• THERMAL-STRUCTURAL OPTIMIZATION 
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TRANSIENT THERMAL ANALYSIS TIME REDUCED BY IMPROVED SOLUTION METHOD 

One Important element of the Langley research program has been the Implementa- 
tion and verification of efficient solution algorithms. A promising set of implicit 
algorithms denoted GEARIB has been installed in SPAR. The desirable feature of the 
6EARIB technique is the ability to adaptively vary the step size throughout the 
temperature history. As indicated in figure 8, various algorithms were used to 
trace the 3500-second temperature history in a section of the Shuttle orbiter wing. 

Two models were considered: a two-dimensional section through the wing depth and a 
one-dimenstional plug through the center of the two-dimensional model. Calculations 
for each model were carried out using explicit Euler and implicit backward differences 
as well as GEARIB. The explicit algorithm is burdened by the need to take small time 
steps (i.e., 0.1 seconds) to avoid numerical instability. The backward differences 
algorithm uses a larger but fixed time step of 1.0 second (determined by accuracy 
considerations). The GEARIB algorithm, by adaptively changing its time step to as 
much as 218 seconds, obtained solutions with significantly less computer time than 
the two previous methods. Additional results (not shewn) were obtained using a 
lumped-parameter thermal analyzer (MITAS) and indicated that for consistent models, 
solution methods, and accuracy levels, SPAR and MITAS solution times are comparable. 
Thus, it is expected that use of GEARIB in this type of analyzer would lead to 
efficiency improvements similar to those obtained in the finite element program. 
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I<educed-basis method for transient thermal analysis 


Along with Implementing Improved solution algorithms such as GEARIB, work has 
been Initiated to develop approximate analysis techniques which have potential for 
significant reductions In solution effort. One such method Is the reduced-basis 
technique which combines the classical Raylelgh-Rltz approximation with contemporary 
finite element methods to retain modeling versatility as the degrees of freedom are 
reduced (ref. 3). The effectiveness of the method depends upon representation of 
local temperatures by a few modes or basis vectors. The reduced-basis technique has 
been successfully applied to the problem Illustrated In figure 9. This problem 
represents a 58-1n. segment of the loafer surface of the Space Shuttle wing and 
consists of a 119-mil. -thick aluminum skin covered by 1 In. of Insulation. The 
combined structure was modeled with two-dimensional, finite elements with 84 node 
points (84 degrees of freedom). A heat pulse reasonably representative of Shuttle 
reentry was applied to the surface and produced temperatures where radiation becomes 
appreciable. Additionally, the thermal properties of the insulation are nonlinear 
functions of temperature and ambient pressure so that the heat transfer equations 
are highly nonlinear. A total of 23 thermal mode shapes were selected from 
solutions of two thermal eigenvalue problems: the first based on material 

properties evaluated for uniform temperatures of SGO^R, and the second based on 
temperatures from a steady-state problem with averaged heating and thermal 
properties. The resulting temperatures are compared with temperatures obtained from 
a SPAR thermal analyzer solution of the full system of 84 equations. The 
temperature Histories shown on the figure agree very well and Indicate that the 
reduced-basis technique can approximate temperatures from the full system within 
20®R over the entire heat pu’-e. Efforts are continuing to find the minimum number 
of basis vectors for acceptable temperatures and to demonstrate the technique for 
larger problems. 


SURFACE HEATING PROFILE 
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ORBITING TRUSS SPACE STRUCTURE FOR TESTING INTEGRATED FINITE ELEMENTS 


.* _ Development of Integrated finite elements has recently been focused on elements 

* with radiation heat transfer (ref. 4). The next three figures reproduced from 

reference 4 describe some recent results from that work. The test problem for 
evaluating the integrated elenents is shovm in figure 10. A three-member module of 
Ir an orbiting space truss is shown. A typical truss number receives solar, Earth- 
- emitted and Earth-reflected heating and emits thermal energy to space. In a 

geosynchronous orbit, solar heating predominates. Incident normal flux to a truss 
member varies significantly as a member changes orientation with respect to the 
vector. As the orbiting structure enters and emerges from the Earth's 
shadow, significant changes in incidert heating occur. Member temperatures and 
structural deformations depend strongly on the time-dependent heating and member 
material and surface properties. 






COMPARATIVE TEMPERATURE DISTRIBUTION OF A THREE-MEMBER ORBITING TRUSS 


Temperature distributions from conventional and integrated finite elements for 5 • 
the three-member truss are shown in figure 11 for a typical orbital position of the 
truss space structure. The term integrated element refers to a thermal element • 

formulation using a nodeless variable to obtain a quadratic temperature variation in 
a two-node element so that accurate thermal forces may be obtained for thermal i 

stress analyses using two-node structural elements. Solutions were obtained for * 

models using a single conventional and a single integrated finite element to ’ 

represent each truss member. A solution using 10 conventional elements per member 
is used as a baseline for comparing conventional finite elements with integrated ? 

finite elements. The conventional elements did not give a good representation of ■ 

interior member temperatures. The nodeless variable elements predicted member • 

temperature distributions and nodal temperatures very accurately with small r 

deviations from the reference solution. i 
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COMPARATIVE DISPLACEMENTS OF A THREE-MEMBER ORBITING TRUSS 

Histories of a typical member elongation for one orbit are compared for the 
various finite element models (fig. 12). The member elongation as computed from the 
one element per member conventional model temperature distributions show up to 44 
percent deviation from the reference solution. These discrepancies arise because 
the linear cenperature distributions predicted by the conventional elements give a 
poor representation of the averao? member temperature. Since the nodeless variable 
elements represent temperature d: ;.tributions very accurately » the elongation 
predicted from these temperatures shows excellent agreement with the reference 
solution with the largest d1scre(.-.ncy less than 1 percent. 
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OF ISOTHERMAL TRUSS ELEMENT 

Although the rather rigorous Integrated elements discussed In the previous 
figure are In the research stage and hold promise for future applications, less 
rigorous state-of-the-art analysis procedures are being applied which employ 
engineering assumptions. For example, several modeling simplifications are often 
utilized In thermal analysis of space trusses (ref. 5). For finite element analysis 
of space trusses, the simplifications are embodied In the Isothermal truss element 
(fig. 13). The Isothermal element neglects conduction, radiation exchange with 
other elements, and Interelement shadowing. The element temperature is determined 
by a balance of radiation heat transfer between solar and Earth-reflected (albedo) 
energy and radiation from the element to space. The figure typifies the temperature 
distributions obtained with the Isothermal element. Temperatures around the 
triangular substructure of the orbiting truss are compared. The exact temperature 
distribution Is given by the piecewise linear representation. The Isothermal 
elements obtain good overall element temperatures but miss the details near the ends 
of the element. Consequently, use of Isothermal elements can be useful for overall 
temperature distributions and deformations but would not be expected to be 
sufficient for detailed deflections or stresses, particularly at the joints between 
trusses. 



ENERGY 
EMITTED 
INCIDENT 
HEATING 

• CHARACTERISTICS 

NEGLECTS CONDUCTION. INTERELEMENT 
RADIATION AND INTERELEMENT SHADOWING 
GOOD AVERAGE MEMBER TEMPERATURES 
PERMITS SIGNIFICANT MODEL SIZE REDUCTION 
MISSES DETAILED DISTRIBUTIONS 




EVALUATION OF APPROXIMATE THERMAL MODELING 
AND ANALYSIS TECHNIQUES FOR SPACE TRUSSES 

As Indicated previously, thermal analyses of space trusses have often relied on 
simplifying assumptions as«:oc1ated with the isothermal element concept. As 
Indicated In figure 14, there Is a need to evaluate the errors associated with a 
number of simplifications Including the neglect of conduction, interelement 
radiation and shadowing as well as the use of quasi -steady analyses (neglect of 
transient effects). In addition, two approximate analysis techniques will be 
evaluated: namely, the reduced-basis method (fig. 9 and ref. 3) and a thermal modal 

decomposition method. The test problem being used for the first set of evaluations 
Is a truss model obtained from reference 5. 


• APPROXIJMATIONS TO BE EVALUATED 
• MODELING APPROXIMATIONS 

• NEGLECT CONDUCTION 

• NEGLECT INTERELEMENT RADIATION 

• NEGUECT INTERELEMENT SHADOWING 

• ISOTHERMAL ELEMENTS 
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LONG DURATION EXPOSURE FACILITY RADIATION THERMAL ANALYSIS 

Along with many of the complicating factors associated with the thermal 
analysis of the Shuttle orbiter (fig. 2), radiation heat transfer effects have a 
large impact on analysis needs for large space structures. Recent experiences at 
Langley In computing view factors for a future Shuttle payload demonstrate the need 
for efficient calculations of radiation view factors. The Long Duration Exposure 
Facility (LOEF) Is depicted in figure 15. The model shown Included 440 panels and 124 
occluding surfaces, and required 60,000 view factors. Thermal analysis of this 
structure required 32 hours of CPU time using an in-house developed computer program 
on a COC 6600 computer and indicates that radiation analysis of large space 
stuctures such as a space station is likely to severely strain current techniques. 

The results of this analysis along with other similar experiences for Shuttle 
components suggest the need to take a look at more efficient programs such as TRASYS 
(ref. 6} and to attempt to improve techniques used for radiation analysis. 
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RADIATION THERMAL ANALYSIS RESEARCH 

Since complex space station geometry and heating-induced geometry changes will 
strongly Influence radiation heating and interelement radiation exchange* Improved 
methods for calculating view factors and heat fluxes are needed. The TRASYS program 
has undergone significant upgrading as described la reference 7 and will be 
evaluated. This program is used extensively at Johnson Space Center for Shuttle and 
spaoicraft applications. Langley Is also developing Improved radiation capability 
under a grant with the University of Washington. As Indicated In figure 16, 
improved numerical and code capabilities will be demonstrated In a study code and 
will Include such effects as shadowing, reflections, reradiation, adiabatic 
surfaces, planes of symmetry, solar flux, and penumbral effects. Techniques 
currently being considered to Improve efficiency are double area Integration, 
contour Integration, and Monte Carlo methods. After the Improved capabilities have 
been demonstrated for small problems, the SPAR thermal analyzer will be updated to 
Include the improvements. The goal Is to develop in SPAR the capability to compute 
radiation view factors and fluxes, then temperatures and deformations, and finally 
feed those deformations back to the view factor calculation to obtain the view 
factors for the deformed structure. Results from part of this effort may also be 
Incorporated In lumped-parameter thermal analyzers, and some products of the research 
could be Incorporated In TRASYS or TRASYS- type programs. 


• OBJECTIVE 

• IMPROVE EFFICIENCY TO CALCULATE VIEW FACTORS AND FLUXES 


• APPROACH 

• DEVELOP IMPROVED NUMERICAL AND CODE CAPABILITIES 

• DEMONSTRATE WITH STUDY CODE 

• INCLUDE SHADOWING. REFLECTIONS. RERADIATION. ADIABATIC 
SURFACES. PLANES OF SYMMETRY. SOIAR FLUX. PENUMBRAL 
EFFECTS 

• TECHNIQUES 

• DOUBLE AREA INTEGRATION 

• CONTOUR INTEGRATION 

• MONTE CARLO 

• RESULTS 

• EFFICIENT NUMERICAL AND COMPUTER TECHNIQUES 

• INCORPORATION IN SPAR THERMAL ANALYZER 

• TRANSFERABLE TO LUMPED-PARAMETER ANALYSES 

' ' ' .Figure 16 
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ROLE OF ADVANCED COMPUTER HARDWARE 

Discussions Involving advanced computer hardware at a recent NASA symposium. 
Computational Aspects of Heat Transfer In Structures (ref. 2), brought out several 
observations on the role of advanced computer hardware (fig. 17). Significant 
reductions in run times were achieved for the unmodified SPAR thermal analyzer on a 
CYBER 203 and for SINDA on the CRAY. The discussions also indicated that the effort 
required to restructure (vectorize) large operational codes to take full advantage of 
vector processing requires such a large personnel Investment that It is not generally 
cost effective. There might be justification for building a vectorized program from 
the beginning provided the community of potential users was significant. Benefits 
mentioned for minicomputers were their Interactive capability and virtual memory which 
permits solution of large problems. Benefits mentioned for the supercomputers were 
their large In-core capacity which permits rapid solution of extremely large problems. 


• SIGN I FI CANT REDUCT I ON IN RUN TIME ACHIEVED (SPAR ON 
CYBER-203. SINDA ON CRAY) 

• RESTRUCTURING PROGRAMS FOR VECTOR COMPUTERS NOT 

COST EFFECTIVE . , . 

• BENEFITS FROM MINIS (INTERACTIVE. VIRTUAL MEMORY) 

• BENEFITS FROM SUPER COMPUTERS (LARGE CAPACITY) 
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THERMAL ANALYSIS OF 30 -METER PRECISION DEPLOYABLE ANTENNA 




One of the early examples used to evaluate finite element thermal analysis on 
the CYBER 203 computer was the antenna structure shown In figure 18. This was a 
fairly small, model , containing only 55 grid points and 183 elements, but did have the 
complexities of Interelement radiation and shadowing. A transient analysis was 
performed for one oeosynchronous orbit and solution times on the CYBER 175 and 
CYBER 203 computers were compared. The 31-percent savings of computer time was 
considered significant but not entirely satisfying since greater savings had been 
anticipated. It was decided that since the program operated In the scalar mode and 
the program had not been restructured to take better advantage of the vector 
processing capability of the CYBER 203. the time saving was as much as could be 
expected. This reasoning, however, led to a subsequent task of performing some 
vectorized programming to assess the benefit of vector processing by the CYBER 203. 




SOLUTION TI/VtTsec 

CYBER 175 

CYBER 203 

46.4 

32.1 


•scalar mode 


• 55 GRID POI^^rS 

• 183 ELEMENTS 

• INTERELEMENT RADIATION 

• TIME- DEPENDENT SHADOWING 

• TEMPERATURE HISTORY FOR 24 HOURS 
OT = .01 hr 
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IMPACT OF VECTORIZATION ON SOLUTION TIME ASSESSED 

A pilot program was written to determine the benefits of the vector processing 
capability of the CYBER 203. The program solves for the transient temperatures in 
an insulated cylinder modeled with SPAR solid thermal elements. Figure 19 shows the 
vectorization stages that have been completed. The CYBER 203 run time in the scalar 
mode is shown in the first line of the table. The next entry displays the benefits 
from obvious conversions of DO loops to explicit vector calls and the 'factorization 
of scaling the element conductivity matrices. The subroutines which factor and 
solve the symmetric banded system of equations were replaced by a vectorized 
subroutine from the CYBER 203 system math library. The answers produced were 
identical, and the tiir<e required for this operation was cut by almost two-thirds. 

The library routine uses a vector length of half the bandwidth plus one. For the 
insulated cylinder model this was 26. A larger bandwidth would produce more savings 
here. The single most time-consuming operation is the multiplication of the 
conductivity matrix, denoted K, by the temperature vector, denoted by T. The 
reprogramming of this operation so that relatively long vectors were operated on 
(ref. 8) was responsible for the largest time savings and brought the total solution 
time down to 33 seconds, a reduction to nearly one-third the original solution time. 


PROBLEM: 1000 sec TEMPERATURE HISTORY OF 800 NODE CYLINDER. DT = 2.0 sec 







LEVEL OF VECTORIZATION 

CPU 

TIME 

ORIGINAL “ NO VECTORIZATION 

92 

EXPLICIT VECTOR CALLS FOR 
OBVIOUS LOOPS AND SCALING 

85 

VECTORIZED ROUTINE FOR 
SOLUTION OF EQUATIONS 
(MATH LIBRARY ROUTINE) 

65 

VECTORIZED K T OPERATION 

33 
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PROPOSED APPLICATIONS OF HEAT PIPES TO AEROSPACE STRUCTURES 


The u$e of heat pipes for reducing temperature gradients and the associated 
deformations and stresses has been proposed for a number of applications (fig. 20) 
such as a beat- pipe-cooled leading edge for the Shuttle Orbiter wing (ref. 9) and 
a heat pipe sandwich panel radiator concept for the orbiter payload bay doors. For 
large space antennas* a heat pipe sandwich panel has been proposed for the reflector 
structure (ref. 10). Finally, heat pipes have been proposed for thermal management 
on space station concepts (ref. 11). These’ proposed applications lead to the question 
of whether existing capability is adequate to handle detailed analyses of heat pipe 
thermal performance. 





CONSTRAINT TRACKING FOR TRANSIENT THERMAL RESPONSE 

An Important aspect of space structure analytical design will be the 
requirement that time- dependent temperatures and stresses be less than specified 
allowable values (constraints) for all appropriate values of time for all orbits. 

The traditional approach for such problems has been to calculate response quantities 
for a few "time slices" and design the structure so that constraints are satisfied 
for each time slice (fig. 21). An extension of this approach called "discrete 
times" calculates response quantities at many discrete times which span the time 
domain (ref. 12). The response for each discrete time Is treated as a separate 
constraint and the structure Is designed to satisfy all the constraints. For the 
example shown, two such constraints are not satisfied and the structure must be 
resized. A third method to handle time-dependent constraints known as the "critical 
times" method Is to calculate times when the largest response for a given design 
occurs, size to satisfy the constraints for those times, and periodically update the 
critical times as the optimization proceeds (ref. 13). 


APPROACH 

• USUAL TIME SLICING 


• DISCRETE TIMES METHOD 

(ref. 12) 


DESCRIPTION 

CALCULATE RESPONSE AT DISCRETE TIMES 
TREAT EACH AS A LOAD CASE 
DESIGN TO SURVIVE EACH LOAD CASE 

" CALCULATE RESPONSE AT DISCRETE TIMES 
TREAT RESPONSE AT EACH TIME AS A CONSTRAINT 
SIZE TO SATISFY ALL CONSTRAINTS 


• PREDICTED CRITICALTIMES 
METHOD 

(ref. 13) 


CALCULATE TIME (S ) WHEN LARGEST RESPONSES 
OCCUR 

SIZE TO SURVIVE AT THOSE TIMES 


RESPONSE 



DESCRETE TIMES METHOD 
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OPTIMUM DE5IGN FOR THfRMAL DISTORTION CONSTRAINTS 

k *k*'*c®® flexible space structures sur..; as antennas and optical systems typified 
by the Space Telescope (fij. 22) require the consideration of extremely tight 
constraints on temperatures, temperature gradients, and thermal deformations. For 
example, part of the optical control system In the Space Telescope is permitted no 
more than a 0.0005 F temperature difference over an 18-inch span. Typical methods 
for controlling temperatures and deformations in the telescope include the use of 
theimal coatings and heaters. A contractual effort is. being initiated with the 
Perkin-Eimer Corporation to adapt mathematical optimization procedures to the 
problem of automatically selecting parameters such as coating locations, heater 
locations, and sizes to meet specified tolerances on the deformations and 
temperatures of orbiting structures. This effort will build on the work of the 
principal investigator's doctoral dissertation (ref. 14) which has already been used 
to influence design procedures for the Space Telescope. As part of this work, 
thermal analysis techniques will be modified as required to assure that the hioh 
precision needed in computed temperatures is achieved. 
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SUMMARY 


This paper has Intended to discuss a number of Issues and needs relative to 
thermal analysis of large space structures and space stations (fig. 23). Some 
Indications of trends In the Langley thermal -structural analysis research program 
consistent with the Issues and needs are also presented. The main heat transfer 
mechanism In space Is radiation; consequently, there Is a need for a strong thrust 
on Improved radiation analysis capability. Also the Important Interactions among 
temperatures, deformations, and controls need to be accounted for. Finite element 
analysis capability seems to be lagging behind lumped -parameter capability for heat 
pipe analysis. Briefly, the Langley plan will Include Improving radiation analysis 
capability, evaluating the errors Involved In certain approximate analysis and 
modeling techniques for large space trusses, and continuing the development of 
Integrated thermal -structural finite elements with an emphasis on radiation heat 
transfer. Work will be Initiated to develop finite element analysis techniques for 
heat pipes. Finally, optimization research activities will be oriented toward 
methods to design flexible orbiting structures to account for thermal and thermal- 
deformation requirements. 


• THERMAL ANALYSIS NEEDS OF LSS PRESENT RESEARCH OPPORTUNITIES: 

• RADIATION ANALYSIS 

• INTERACTIONS: TEMPS. DEFORMATION. CONTROLS. ETC. 

• HEAT^PI PE MODELING. ANALYSIS. OPTIMIZATION 

• PLANS 


• ORI ENT THERMAL-STRUCTUAL ANALYSIS ACTIVITY TO LSS 

• MAJOR THRUST ON RADIATION ANALYSIS IN SPAR 

• EVALUATE SIMPLIFIED MODELING AND ANALYSIS TECHNIQUES 

• CONTINUE INTEGRATED FINITE ELEMENT DEVELOPMENT 

• INITIATE FE HEAT-PIPE MODELING AND ANALYS I S RESEARCH 

• INITIATE OPTIMIZATION APPLICATIONS TO LSS 
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la feedback control problems It has recently been established that there Is 
an optirmm contTolter order for a given mission (set of control objectives) [1]. 

The conclusion is that it is possible to make the model too simple (something every- 
one knows), and it is possible to nuke the model too oomplexi The engineering expla- 
nation for this phenomenon is that models of too low order have large effects from 
errors in model order, whereas models of too high order have large effects from 
errors In parameters. One might chink Chat parameter identification methods could 
be added to remove this latter difficulty. But, alas, it is well known in identifica- 
tion theory [2]-[9] Chat overparameterization actually degrades the performance 
(convergence) of the Identification algorithms. Hence, further research is needed 
to determine methods to find the optimum controller complexity (controller order and 
y. controller parameterization). That is, Che modeling problem and the control problem 
must be united within some Integrated systematic approach. 

The method currently under Investigation for accomplishing this is component 
cost analysis [1]. Component cost analysis determines Che critical components in a 
dynamic system, and varied applications of these concepts lead to an algorithm which 
integrates the following design problems. 

• Which components of the system structure should be redesigned and what 
parameters of the redesign are Important to change? 

• Which sensors and actuators should be redesigned; where should the 
sensors and actuators be located; what type of sensors and actuators 
should be used? 

• Which parameters will be most critical to identify? 

• What controller order and dssociated -optimum controller banduidth is optimum? 

• What sample rate should be selected for digital control? 


1.0 INTRODUCTION 



Flexible structures and their dynamics have been studied for over a century. 
However, only recently has there been an interest in the active control of flexible 
structures. Such Interest was piqued in the 1960’s by a flexibility-induced insta- 
bility in the USA’s first satellite [10] and, more recently, by sophisticated require- 
ments for precision-controlled structures in space for astronomy, communication 
networks, near-Earth scientific studies, and space solar-power alternatives [11]. 

The rapid development of computers and control theory in the 1960 's has encouraged 
active control applications for other structures as well, such as flutter suppression 
in aircraft [12], and active damping of bridges and tall buildings [13]. This is 
not to say that active control is needed in every structure, however, and there is 
no clear way to make the decision of when and how much control effort is needed in 
a structure. There is a need to study the dynamical properties of the mechanical 
system with a view toward discerning what Improvements in performance can easily be 
made by redesigning the structure and what improvements must be left for active con- 
trol functions. This beneficial Interaction of the dynamics and control disciplines 
in the development of a rotational design methodology has not yet occurred to any 
mature degree. Usually the structure designs and the control designs occur sequen- 
tially. This luxury cannot be afforded in the future. Stringent requirements 
force us to provide better coordination between structure design, control design, 
and controller software design. 

Some of the reasons that the control of flexible spacecraft can be a difficult 
task are briefly described by the following three problems. 





(i) The Model Error Problem ' 

The space structure Is usually constructed of lightweight materials, and thus 
the assembled structure Is very lightly damped. This uniqueness of light damping 
for the space structure makes the control design extremely sensitive to modeling 
errors, since the slightest perturbation of truncated modes by control action can 
shift these eigenvalues into the right half-plane. Also there Is the usual uncer- 
tainty in the computation of the modal data. This problem is especially critical 
for spacecraft since modal data uncertainties cannot be removed before flight, due 
to the difficulty of testing the extremely lightweight structure In a 1-g environ- 
nent. 

(11) The Limited Controller Software Problem 

The practical limitations of memory and speed of onboard computers mean that 
only controllers of constrained dimension can be considered. These constraints can 
severely reduce the performance capabilities of the controlled system due again to 
the effect of modeling errors Imposed by the controller order constraints. (An 
infinite dimensional system controlled by finite controllers Immediately suggests 
that the standard "optimal" state feedback solutions are not going to be realized) . 
Thus, limited software serves only to compound the model error problem by constrain- 
ing the order of the controller and by adding delays In the feedback loop . 

(ill) The Performance Requirement Problem 

Of course the model error problem and the limitations of software pose no seri- 
ous threat to the mission If the performance requirements are quite lenient. Thus, 

. the degree to which (1) and (11) pose problems is directly related to the severity 
of the performance requirements. Therefore, early research on the subject sought 
to help with the trade-offs between performance and modeling errors (including those 
Induced by controller software limitations). 

2.0 MODELS OF SPACE STRUCTURES 


Those portions of the structure resembling beams, plates, and membranes might 
reasonably be idealized as a material continuum. The resulting partial differential 
equations (PDEs) contain all the modal data over an infinite spectrum [l4]. Other 
parts of the structure might contain trusses or complicated connections which re- 
quire a finite element formulation of the model, resulting in a set of ordinary 
differential equations (OOEs) [l5]-[l7]. Also, the dynamics of actuators and sen- 
sors are usually described by ODEs. This combination of distributed- parameter models 
(PDEs) and lumped -parameter models (ODEs) must eventually be reduced to a finite 
set of ODEs. The discretization of the PDEs iiiust be accomplished so that the fre- 
quency spectrum over which each of the subsystems (actuators, sensors, substructure 
1, substructure 2, etc.) are modeled Is consistent. Otherwise, troublesome dynamical 
interactions between subsystems might be unintentionally concealed at the outset. 
Thus, each of the substructure models might be truncated prior to the assembly of 
the composite model (2.1). 


I structures, (2.1a) 
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where the Inertia or "wass" ^trlx H' ■ M’ >0 is positive, definite, and synnnetrlc; 
the stiffness matrix K* ” K'* ^ 0 Is positive, semidefinite, and symmetric; the inter- 
nal energy dissipation due to damping Is q^C*q and C » C'^ i 0 is positive, semi- 
definite, and symmetric; and any gyroscopic term is due to 6’ •• which is skew- 

symmetric. The disturbance models assumed for W., V., w , V are all zero-mean white 

“ “ * S a N o'"a 

noise processes. The dimensions of the vector are X^e R,Xj^ER,qeR,U^EK“, 
U e R***, y E R y, Z E R*^Z. xhe coordinate transformation 


q = Tn 

is often made to put the structural subsystem coordinates 


(2.2a) 


n + Z,n + Z2n = 8U3 . 8 = T^8*. Zg = . Z, = T^[C + G']T (2.2b) 
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P » P'T . R = R*T 


(2.2c) 


Let the first n^ elements of q and n be associated with rigid-body (zero frequency) 
modes of the structure. Then these matrix partitions apply 

Zn = t"^K'T = block dlag [0 = dlag [w?, ..., w?] = Tg ^'Tg (2. 2d) 


T^M'T = I , T = [T^ Tg] , T'CT = block dlag [0 a] , A ^ 0 
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ij = CnJ nj] , R = [Rr 


= [Pr 


(2.2e) 


(2.2f) 


«rhere the zero in (2.2e) holds if the rigid modes are undamped, and 6^^ = t 

• Now (2.2) becomes 
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3.0 ON THE STRUCTURE OF MODELING ERRORS 
Let the composite model (2.1) be written In state form 


X»Ax + Bu + Dw+f, 
z = Mx + V + g. 


“ (n^. n. xj. xT) . 

, xcR , n i 2n.+a+s (3.1) 

“ ' “I- 


where nonllnearltles f and g might be added to the model for performance evaluation 
purposes, but might be ignored during control design. The associated parameters 
(A,B,0,M) may be deduced from (2.1), (2.2). There may be several stages of model 
simplification between the most general model used for simulation and performance 
evaluation before flight and the model upon which the control design is based. To 
simplify the discussion, we discuss only two models. In this section, model (3.1) will 
represent the physical system (admittedly In this case, the state x Is infinite- 
dimensional and the parameters A,B,D,M are not precisely known). The reduced model 
used for controller design Is 


Vr+ Vr 


(3.2) 


^R “r^'R ''r 


and we postpone the discussion of how (3.2) might be derived. Our current Interest 
concerns the characterization of the differences between any two models (3.2) and 
(3.1). In order to match Identically the measurements z(t) actually obtained from 
the physical system, one could define vector functions of time e (t), e (t) as those 
which drive the state and measurement equations in such a way ^ ^ 


*R ' Vll ♦ V * Vr * 


(3.3) 
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that z(t) matches the actual measurements. 
[lS3 to be composed of four parts 



Such model error vectors have been shown 

(3. A) 


where is due only to parameter errors, is due to errors In model order, 6^ is 
due to neglected disturbances and is due to neglected nonllnearitics. Of course, 
none of these e^f can be' known a prion. Parameter adaptive control 

and identification methods L19il~C2l} strive to drive 6^(t) to zero. Such methods 
can be effective uhen ia the dominant source of error in (3. A). Terms contribut- 
ing to depict Che coupling between retained and truncated equations in the 
(Inflhlce-dlmenslonal) model underlying the physical system. Much attention C223- , 

[24] has been devoted to the reduction of the "spillover" terms e. (and their 
corresponding closed-loop consequences). Such approaches can be effective if e^. hap- 
pens to be the dominating term in (3. A). Orthogonal filters [25], [26] make fewer 
assumptions about e(t), save that it is square integrable. This more general charac- 
terization of model error has the potential advantage of simultaneously accommodating 
errors of the type (3. A), but the method also has disadvantages in the design stage. 
The method requires the selection and storage of a set of Independent functions (to 
be used for fitting the actual error function), and the choice of these functions 

is not unique. One choice is to use elements of the state transition matrix for a 
higher order model as in [27]. Another choice Is to use orthogonal functions as in 

[25] . 

Having described the nature of the modeling errors, we wish now to be more 
explicit about thetr^ effects. These effects wil) now be described, first in terms 
of stability and then in terms of a quadratic performance metric. 


4.0 STABILITY AND PERFORMANCE IN THE PRESENCE OF MODELING ERRORS 


4.1 Stability and Modeling Errors 

Let 

u • G 

U(S) - [G(sI-A|,-B|,G+FM|,)'’f]z($) «•« 

Xr = ArXr + BrU + F(z - MrXr) 

describe the dynamical controller used to control the system (3.1). G is usually 
chosen to stabilize [Ar + Br6] and F is usually chosen to stabilize [A - FM], 
although the controller poles, ^-jCAr + B.>G - Ri], j = 1, ..., r, should also be 
stable [28]. There is no unique relationship between the parameters of the control- 
ler (4.1) and the evaluation model of the system (3.1). The controller parameters 
(An, Bp, Mp) may or may not be related to some reduced-order model obtained from 
(3.1). For example, the simpler model (Ar, Br, Mp) might have been derived directly 
from physical laws, but using an idealization of the system that was simpler than 
the idealization that led to (3.1). For these reasons the phrase "parameter errors" 
has no unique meaning. In fact, none of the terms in the decomposition (3.4) is 
unique. 

*By idealization we mean the set of hypotheses within which the dynamical system is 
assumed to move. For example, different idealizations of the same system might in- 
clude: an elastic material continuum, a set of connected rigid bodies, a rigid body. 



It has been established £293i [30} that 

(a) • When the controller gain GF in (4.1) is "small enough" , the modeling 
errors will not destabilize the closed-loop system which is open-loop 
asymptotically stable. 

(b) • WJien the controller gain GF in (4.1) is " large enough''^ the modeling 
errors will always destablize the closed-loop system. 

(c) • In the absence of structural damping (C = 0 in (2.1a)), a stablizing 
linear feedback controller (4.1) does not exist. 

The Interesting conclusions concerning (c) are that: (1) the undamped structure 
is used most often in model development foi. *'ontrol design, and (11) the small and 
highly uncertain structural damping chat every real space structure will have is 
vhat makes control possible. For successful application, every control scheme ralles 
upon and presumes the existence of a finite amount of structural damping in the 
open-loop system. 

Rate feedback "guarantees'* to add some damping in the system under the presump- 
tion of infinite bandwidth sensors and actuators [30}, but too much rate feedback 
actually degrades output performance (such as the llne-of-slght errors in a tele- 
scope). See [ 31 } for a spacecraft example of this phenomenon. 

4.2 Performance and Modeling Errors 

We cite here the performance available with and without modeling errors. If 
one chooses G and F in (4.1) to be optimal for the model (A^, li|^, M^) with noise 
intensities ^ and Vn respectively for the zero-mean vdilte noise processes ^{^(^) ^<>4 
V|j(t), and £j there Were no model errors e(t) = 0, then the trade-off between fhe 
optimal output performance 

y™ 5 lim E||yj|2 (Q > 0) (4.2a) 

t-*« 

and the optimal input performance 

= lim Eljull^ , R = R o > 0 (4.2b) 

t-*«> o 

(where p is an arbitrary weighting scalar chosen by the designer) would be 
"hyperbolic" In the sense of Fig. 1. 
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The output of the reduced model is yp -- CpXp. Under the assumption!? of contrpllabil 
Ity of (Aj^, B^) and observability of noth*' (^R* ^^R^ [32] 

6 = 1 bJk, 0 = KAr + aJk + cJqCp - gJrG^ (4,3a) 

F = PmJv^^ , 0 = PaJ + A^P + D^W^dJ - FVj^f’’' (4.3b) 

and small p leads to large control gains. In the presence of inevitable modeling 
errors. Fig. 1 Illustrates the eventual (and perhaps rapid) degradation in output 
performance with Increasing control authority (decreasing p) . This result is not 
predicted, of course, by the standard linear quadratic Gaussian (LQG) theory since 
the theory relies upon the absolute fidelity of the mathematical model. Methods 
to combat such model error effects are the subjects of many papers, but satisfactory 
solutions have yet to emerge. For the lower order controllers, the low control 
efforts in Fig. 1 give large errors and the large controls give large errors. 

Thus, there exists an optimum level of control effort which is a function of controller 
order. 

For any performance metric of the form (4.2), the relationship between stabil- 
ity and the value of the performance metric V is established as follows. Denote the 
closed-loop plant matrix by A, and the output matrix by C, and ignore nonlinearities 
f and g. -Then for the oZosecf— loop system 

, X = Ax + t?w . ' ° -0 c (4.4a) 

y = Cx, V = liin Ell^ll^, a > 0, / = (y'*', u^) (4.4o) 

If the matrix pair (A, C) is observable and (A, P) is dlsturbable then V is finite 
only if A is stable [32], Thus, relative stability information is contained in the 
magnitude of U, and we have the result: 

Proposition 1: Stability margins with quadratic performance metrics 

If observability and disturbability are properties of the closed-loop system (4.4) 
then ^HS serves both as a stability margin and a quadratic perfomance metric. 

Thus, the calculation 

= tr pd^QC , 0 = PA^ + AP + (4.5) 

provides an acceptable performance evaluation if the triple (A, C, P) is dlsturbable, 
observable , There are three problems which prevent the use of proposition 1 in 
guaranteeing stability of physical systems. 

(1) The "physical" system is infinite dimensional. 

(2) Observability and disturbability tests are impossible to do with 
precision on a digital computer, even for a finite-dimensional model. 

(3) The physical system is never observable in the following sense. 
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Proposition 2: Unobservability of physical systems 

As mathematical models increase in complexity, describing more and more_ com- 
pletely the dynamical details of the physical process, the model eventually includes 
unobservable states with respect to any finite dimensional output vector y. 

Prorosltlon 2 explains the sense in which physical systems are never completely 
observable. This proposition is intuitively verified by imagining that microscopic 
phenomena such as molecular motions are not going to be observable in the rate gyro 
measurements. Even though one may not be interested in such trivial examples as 
molecular motions, the useful point of propositions 1 and 2 is that as far as 
stability is concerned there is no clear way to know when minute motions become 
"nontrivial'*. Thus, the propositions serve to remind us that one can never guarantee 
stability of the physical system. This point may be" made from another view. Since 
stability is a mathematical property of a mathematical model. Interpretations of 
stability for the underlying physical system must be accompanied by precise state- 
ments of the type of modeling errors ignored. Thus, the term "stability guarantees" 
refers only to properties of the assumed model, and these comments must be kept in 
mind when reading the "stability" results of this (or any other) report. 

5.0 COMPONENT COST ANALYSIS 

The central Idea to be discussed in this section is to exploit the precise 
statement of the optimal control p; oblem in order to predict which system components 
will make the largest contributions in the total quadratic performance criterion. 

These components are retained, and the balance is discarded from the model. A 
"component" of tl.u system can be any subset of coordinates designated by the analyst. 
For example, all of those coordinates associated with a particular substructure 
(antenna, solar panel, etc.) might be called one "component" for the purpose of 
assigning a value to the component for its contribution ip the total performance 
criterion. The general ideas of such "component cost analysis" (CCA) are described 
in [l, 33-36]. As another example, each modal coordinate might be designated as a 
"component". In this event the procedure becomes modal cost analysis (MCA), and the 
result is that the contribution of each ">ode is determined for the given quadratic 
performance criterion. MCA has been applied to models of flexible structures in [37]. 

5.1 Perfcrmance Objectives 

The performance of the dynamical system must be judged by a specific criterion. 
One may require specific pole locations, focusing only on stability, but that leaves 
the question of eigenvectors and output performance imprecisely specified. Stabil- 
ity is clearly not a sufficient design goal. The linear quadratic Gaussian (LQG) 
problem has the advantage that poles and zeros are Involved in the design specifi- 
cation and that the motion of specific variables of interest can be penalized direct- 
ly by inserting these variables into the cost function to be minimized. On the other 
hand, these specified variables may be regulated satisfactorily while other variables 
become unbounded. The earlier Proposition 1 states that observability is sufficient 
to prevent this situation, but it may not be known whether all the potentially unstab- 
le motions are observable in the selected variables for minimization. Thus, in 
oversimplified terms, the "classical" approach is to "design for stability" (and then 
we must check for performance) , whereas the LQG approach is to "design for perfor- 
mance" (and then we must check for stability) . Insights into the best of both 
methods are required for successful designs. In this section the LQG methods are 
presented. 
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We define by the vector y the collection of all variables we wish to directly 
control. For example, if attitude control of the^rigid body [described by n in 
(2.3a)] is the objective, then, one might chgose y = rip. Alternately, if vibration 
suppression is the only objective, then y = ng wight be chosen. The choice of a 
weighting matrix Q in (4.2) ^s often dictated by energy considerations. As an 
illustration the choices {Q = 1, y - lead to y'Qy - kinetic energy in elastic 2 
moges. The potential energy in elastic inodes is characterized by the choice {Q ^ u , 
y = Hg}. The expected value operator E is required in (4.2) due to the presence of 
random "noisy" disturbances in the actuators and sensors (2.1). The control mean- 
square effort ujjs is added to the performance metric, where the scalar p is an 
arbitrary weighting scalar which trades off control performance y^g versus control 
effort U|^S* practice, p must be chosen so that U{u]5 does not exceed the physical 
bounds of the actuators or structure load design constraints. 

The objective of component cost analysis is to decompose the performance metric 
V.lnto "component costs" I/., 

= I I/,- (5.1) 

1 = 1 ^ 


V = 


^MS 


P U| 


MS 


I \ 
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where Oq is the number of components in the system and is in the in situ contri- 
bution of component i. We choose in the next section to define "components" of the 
epen— loop system as "modal" coordinates and their sensitivities. In section 5.3 
we choose to define a "component" of the closed-loop system as a state of the 
dynamical controller. In each case we intend to assign a relative importance to 
each component by determining its ranking in the manner 


‘'l i ‘'3 - ••• - 


(5.2) 




and truncate components with small component costs V- from the system. The next two 
sections give the necessary mathematics. 

5.2 Modal Cost Analysis 


Ignoring the sensor and actuator dynamics in this section, and considering the 
nongyroscopic system. 




Hg + APg + m Pg = BgWg, 


Pg t R 


N 


EWg(t) = 0, EWg(t)Wg(t) = Wg6(t-T) 
V = y = lim E||yl|2 

t-K«> ^ 


(5.3a) 

(5.3b) 

(5.3c) 






where only the noisy part of the control forces lljj = U + Wg is considered in this 
section. When N is small enough for closed-loop control calculations, we will not 
Ignore the effects of U in model reduction decisions. However, we now presume that 
N is very large in (5.3a) and U cannot yet be computed. 

Our Immediate objective is to ascertain the contribution of each mode of the 
system in the overall cost value (/. We presume at this point to be dealing with a 
system of very high order whose control inputs u(t) are not yet known. Yet to 
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Ignore the source of excitation altogether would be a mistake since the final judg- 
ment of the quality of the model is in the presence of actuator activity. Until more 
is known about u(t), (after control design considerations) we assume that U = U„(t) 
is a white noise of random disturbances being propagated through the actuators due to 
electronic noise in the electrical or magnetic amplifier devices, etc. The consider- 
ation of the control and its effect in the model reduction process will be postponed 
to the next section. To proceed we now need the following definitions. 

Definition: Let be any subset of state variables of the linear system 

X = Ax + Dw , Ew(t) = 0 , Ew(t) w^(t) = WS(t-T) 
y = Cx 


= (... x.^ .. . ) , X. , xeR"’' , = y n. 

1 1 x 1 


(5.4) 


4 AT 

1/ = lim E Y , Y = y‘0y 


Here x^ will be called the ith "component-*' of the system and the "component cost" 
associated with component is defined by 


A 1 aV 


3Y A ^ 3Y , 

3x. " 'ax. * •*•* ax. ‘ 


(5.5) 
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in. 


It is easy to verify that the total cost value V is the sum of all component costs 


i 


\i = 



and that the component costs may be computed by 


(5.6) 


= tr[C^QCX]^^ , 0 = XA^ + AX + DWD^ (5.7) 

where [c"^QCX]-j-j denotes the n.xn. matrix partition of c'^^QCX. The above analysis 
(5. 6)- (5. 7) of a linear system is called component cost analysis (CCA) [1). The 
"components" of a system might be defined from physical or mathematical considera- 
tions. From physical considerations x^ might represent any physical component of 
the system such as the states associated with: 

(a) A substructure of the flexible spacecraft (an antenna, a solar panel, a rigid 
body, etc.) 

(b) An electrical or electromechanical element in the system (actuators, sensors, 
amplifiers, filters, etc.) 

From mathematical considerations x^ might represent any mathematical component of 
the state in any transformed (.nonphysical ) coordinates. One such example which is 
common is to define the "components" to be "modes" of the systems. In this case the 
"component costs" (5.7) are called "modal costs" ll], [36], [37]. We choose to 
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now examine the system using such modal cost analysis (MCA)* 

Theorem [30] 

For the system 

•• • ^ ■ * 

+ u. = 0 , {n^(O), n^(0)} specified, i = 1, N 

N 

y = 1 (p^n^ + 

i=V ’ ^ ^ 

the cost funotion 

.. A 


= r y^Qydt = X V. 
Jo i=l ’ 


decomposes into the modal costs given by 

A 1 




'i ' r ^ 4:7 niPi'i q * 

= (n^. n^) 


in the limit as 0. 

Theorem [30] 

For the system 


Hi + 2c;(ij. n- + = b. u J + dJw , W'^N(0, W), u N(0, U) 

® ililllWI w 

N 


y = (Pjn, * 


(5.8) 


(5.9) 


(5.10) 


(5.11) 


(5.12) 


where and w are zero mean uncorrelated white noises with intensities U and W, 
the cost fwiation 

V = Tim E{y‘0y} = 7 V. 

t-«. i=l ^ 

decomposes into the modal costs given by 

- i r- r^'L 


{|^x4 = UlPillq + “^^l|rj|][|lb^llu + lld^llw] ^ (5.13) 


\ = /Oy 


(5.13) 


The proof of the first theorem may be found In [7] and the proof of the second may be 
found In [8J. The value of (5-13) is that both disturbance and control points of 
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excitation are considered in the_subsequent model reduction decisions. Of course, 
the real purpose for including in (5.11) is to anticipate some sort of excita- 

tion through the actuators prior to the actual design of the control law. 

According to modal cost (5.13) the importance of a mode is determined basically 
by the product of three properties of the mode, 

- (time constant) (observabiiity)(disturbabnity + control) ability) (5.14) 

where mode 1 is: unobservable in position if and only if ||p^|| = 0, unobservable in 

rate outputs if and only if ||r^. || = 0, undlsturbable from w if and only if ||d^. || = 0, 
and uncontrollable from if and only if U b. || = 0. The expressions (5.10) and 
(5.13) are general modal costs for matrix second-order systems which have no gyro- 
scopic terms. 

5.3 Component Cost Analysis as a Spacecraft Control Design Algorithm 

The results are now combined to form a design algorithm for flexible spacecraft 
control. Each step in the algorithm is discussed and motivated in some detail. 


The system under consideration has the form (2.1) with 6' 0 and (2.6) holds- 

Kate sensors are collocated with the m actuators. An actuator applies a force 
between two points in the structure, and a sensor measures the time rate of 
change of the resulting rectilinear displacement. Alternatively, an actuator 
applies a torque between two points., in the structure, and. a sensor measures the 
.time rate of change of the resulting angular displacement. 

A number of the m actuators equal to the number of rigid body modes are located 
so that the rigid body modes are controllable. 

There are n inertlally-referenced position sensors located so that rigid modes 
are observaole. 

The number of elastic modes calculated for the structure Is N. 

The largest Liapunov equation of the form (5.7) which can be solved reasonably 
accurately on the available off-line computer has dimension 2N|^ x2N|^. 

The largest Rlccatl equation of the types (4.13), (4.14) which can be solved 
reasonably accuiately on the available off-line computer has dimension Nj^xN|^. 

The largest controller (4.1) which can be accommodated in the on-line computer 
has dimension Nq, where Nq £ Nj^ £ 2Nj^ 2N. 

The highest bandwidth of the available actuators is upj,. 

T 2 2 

The parameters considered uncertain are given by p‘ = (c,, ..., Um| 

R c \ I N • n 

We have in mind the performance measure 


V = lim E{ay^Qy + By^q; + pu\u + \ R 

where for v “ 2N + N m, N_ = i Nn - n_ 

6 6 6 & n P 

To Illustrate the design procedure we use the simple example 


(5.15) 
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\ = B^u + P^w , n^cR 


Hg + 2ciDrig + tD^Hg = BgU + PgW , rtgeR^ 


(5.16) 


y = + PgHg . 2 = M + V • = M^gHg + 


where 


C - 0.005 , u = diag [1,2,5,10], Mp^ = 10, B^ = 1 , P^ = 10, = 1 , Q = 100 

3 = 1, Eu(t)u^(x) = U5(t-x) , U = 1, Evi{t)vi^{r) = W6(t-x), W = 100 
EVp(t)v'^(x) = Vp6(t-x) . 

V„ = 10, Ev^(t)v‘I(T) = V 5(t-x) , = 0.1, = (0.1, 1.0, 0.01, 1.0) 


Bg = (0.1, 0.01, 1.0, 0.1) 

Pg = (0.01, 0.1, 0.01, 0.001) 

and we Cake, for illustrative purposes only, 

N = 4, = 1, m = 1. N, = 3, N„ = 6, N = 2 

i o Lol\ u- Q o 

The measurements Z and Z^ respectively represent position of the rigid body and a 
rate measurement. ^It may be readily verified that all the assumptions (Al-All) are 
satisfied . 

5.4 The CCA Design Algorithm 

Step l ! The Preliminary Model Reduction ; 

Set B = p = Oj = 0 in (5.15) and use MCA (5.13) to reduce the number of vehicle 
modes from (N + n^) to (Nj^ + D^), where is the number of rigid body modes surviv- 
ing the MCA truncation. 

Purpose of STEP I: Reduce the number of modes to a tractable number, but do so with 

knaiiledge of y^Qy, the primary aontrot objective. 

Example of STEP 1: For (5.16), using (5.13) we have 

VCn^) = l/(n^) = 1.00, i/(ri 2 ) = 625.00. V(n^) = 40.40, 1 /( 04 ) " 

Hence, by the modal cost rule for truncation (5.2), mode 4 is truncated and the 
retained modes, listed in order of their modal cost (their predicted Importance in 
the problem), are (Oj., H 2 * hg, n-j). The reduced model is 

»i«.*' 8 u+Pw (5.18a) 


hr “ 8^u + Pj,w 

iig + 2cu'rig + w'^hg = BgU + PgW 


(5.18b) 
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> = PrV ^ * “pr\ * “p • "r ■ “>6 + “r 

cont'd. 

where 

(ij* = dlag [2. 5, 1] , vj = (0.10, 0.01, 0.01) 

B'J = (0.01, 1.00. 0.10) , p; = (1.00, 0.01, o.io) 

STEP II . Rate Feed b ack Design: 

Set u = Un + u„ where 
Or 

where B is chosen- large enough so that Pg"^QPgB-2C(u' > 0 . 

Purpose of STEP II: It follows from theorem 6 In [30] that the control (5.19) is the 

control which is both the optimal measurement feedback control and the optimal 
state feedback, control for 

V = Urn E{f,'^ [2^u.‘ + + uJg/"' u^} (5.20) 

t-K» 

subject to (5.18b). Furthermore, such a control Increases the relative stability of 
all controllable modes and, of course, does not move others. Hence, the system is 
stable in the presence ”of almost all" modal data uncertainties. However, this pro- 
mise is only valid for those modes within the actuator sensor bandwidths. The main 
purpose of STEP II is to increase the damping of those modes that MCA has identified 
as critical to the cost function, and to do this for a larger number of modes than 
the subsequent outer control loop can be optimized for. This allows a control 
spillover "cushion" in the sense that the control spillover from Uq will have to 
push those residual poles (truncated in the Uq design but present in the Up design 
of STEP II) further to the right to destabilize them. Now, in order to be sure 
that STEP II has provided "spillover protection" for the same modes for which such 
protection will be needed later, it is Important that the design of and u be 
"coordinated" to the extent that they are both concerned with the regulation*^of the 
output vector y. The rate feedback design minimizes rates, so choosing Up to mini- 
mize y is the rational thing to do. Hence, we may interpret STEP II in the spirit 
of (5.15) by setting a = a. =0 with U -► U . To match this objective as closely as 
possible with (5.20), we set 

pRq = = 2cm' + (5.21) 

and find the which comes closest (in a least squares sense) to satisfying the 
second equation in (5.21). The result is (5.19). 

Example of STEP II: From (5.18), (5.19) 

G^ = 7.82 (5.22) 

The system (5.18) is now described by 

= B^Uq DpW (5.23a) 
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(5.23b) 


t C2W' * b'g^b;T] + p;» 

y - P^V ' “prV * ''p • V ■ “re"e * V 

STEP III ; Second Stage of Model Reduction ; 

Put (5.23b) in the state form 

X = A„x + B-U- + 0 w 
0 0 0 0 


(5.24) 


y = V 


z = M„x + V 

0 


= {zj . zj) 


where the 2N^ components of X are T = h - - • , Nj_ arranged in any order. 

Compute '* 1 


V. = [c^qc X , 0 = X + AX + D wd]" 
1 ^ 0 o n 0 0 0 0 0 0 


(5.25) 


and delete from (5.23b) those modes with the smaller component costs defined by 


VK . ni ) = l/(n; ) + ^(n; ) 


(5.26) 


where and t^Crig ) are computed from (5.25). The number of equations retained 


■ej 


1 . 


in (5.23b) is - np. 

Purpose of STEP Illr STEP III must reduce the model (5.23) to "Riccati-solvable" 
dimension Nj^ to prepare for design of the control U^. 

Example of STEP III: From (5. 23b)- (5. 26) we have 

‘'(ni , n: ) = 0.60 . l/(n: . ) = 609.20 . l/(n_ . nL ) = -0.02 (5.27) 


1 


-1 


indicating that the equation for Hg (t) is to be truncatec from (5.23b). The model 
now is put into the form 3 


X = Ax + Bu^ + Dw 

0 

y = Cx 

z = Mx + V 


(5.28) 


,, vT - / T ,T ,T -T .,T -,T, 

where X = (n^. Hg^> Hg^ » n^. Hg^. 

STEP IV : Design of the Outer Loop U ^; 

Set = 0 = 8, pRq = G”^ in (5.15) to get the optimal controller 
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(5.29) 


Uq = Gx , X = Ax + BUq + F(z-Mx) 

G = -G^b'^^K , 0 = KA + A^K - KBG^B^K + C^QCa 
F = X[1^V“^ , 0 = XA^ + AX - Xm''’V^MX + DWD^ 

Purpose of STEP JV: The value of R = has been established earlier. Now the 

optimal controller for (5.28) has been computed in (5.29), assuming certain para- 
meters, « 0. Those modes of the structure that are observable in y are also 
observable in It is in this sense that the rate feedback controller of STEP II 

(which focused on the term) is "coordinated". 

Example of STEP IV: The essential data from (5.29) is 

G = (-3.16E-01. -4.41E-02, -5.77E-01, -8.T6E-01, -4.41E-03, 2.62E-01) 

C = (l.OOE+OO. l.OOE+00, l.OOE-01, 0, 0, 0) , a = 7.82 E-03* 


ri.OOE+01 

0 0 

0 0 


0 




Lo 

0 0 

0 l.OOE-02 

l.OOE- 

01 j 



To 0 

0 1 

0 

0 

- 




0 0 

0 0 

1 

0 



■■ 


o 

0 

o 

0 0 

0 

1 





0 0 

0 0 

-7.82E-02 

-7. 

82E-01 




o 

0 0 

-2.08E-02 

-7. 

82E-03 




Lo 0 

-1 0 

-7.82E-03 

-8. 

82E-02_ 










r- 


0 

0 

0 


0 



2.51E+00 

6.35E-02 

0 

0 

0 


0 



2.51E-02 

6.30E-04 

0 

0 

0 


0 


,F = 

2.60E-03 

6.78E-05 

1 

1 . OE+01 

-7.82E+00 


l.OOE+OO 


3.16E+01 

1.62E+00 

l.OE-02 

l.OE-01 

-7.82E-02 


l.OOE- 

02 


3.12E-01 

1 .66E-02 

J.OE-01 

l.OE-02 

-7.82E-01_ 


J.OOE- 

oi_ 


_3.23E-02 

3.36E-01_ 


*Multiplying (5.15) by 1/a yields a parameter p/a to be selected to achieve an accept- 
able control effort. The best value of a is therefore determined from Fig. 1 as that 
choice which corresponds to the lowest point on the dotted curve of Fig. 1. This 
point is found numerically, first picking an a then truncating the controller to the 
desired order Nc, and repeating this process for a new a. Thus, the above value of 
a was found after several passes through STEP V. 
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STEP V: Controller Reduction: 


The optimal system of STEP IV is 

BG 

FM A+BG 


© Fa* bg ~| /^\ + 

Lfm a+bg-fmJ Vx/ |_0 fJ \v/ 

■O'C :](D'“ 


= Ax + Ptv 


(5.30) 


which is now to be evaluated by criterion (5.15). First we show the simpler case of 
known parameters (a^ = 0). 

STEP V-A ; Controller Reduction With Certain Parameters (a^. = 0): 

Delete from the controller 

Nr 


X = [A+BG-FM]x+ Fz , X E R ^ 


u = Gx 

0 


(5.31) 


those states (n^, ) with the smallest costs t'(n^» n^- ) = defined by 

V. = Cc^acx]jj , 0 = Xa"^ + AX + UtJV^ , 0 = block diag [Q, G"^]. ,, (5.32) 

for j = i i = 1, Nj^ vnere V. = V(n^-) and = l/(n^) if X has the form 

X* = (■•• n^. > * ...). This yields the reduced controller 


V ’ “o ^R^R ’ ’'r 


(5.33) 


where Ar is obtained by deleting the indicated rows and columns of [A+BG-FM], Fpj is 
obtained by deleting the same rows of F, and Gj^ is obtained by deleting the same 
columns of G. 

Purpose of STEP V-A: Reduce the order of the controller (5.31) to the order Nq accep- 
table by on-line software limitations. 

Example of STEP V-A; The calculation (5.32) reveals that the cost-ordered components 
of the controller are 


:T 


~ (Op* Dpt 


"eo’ "e, 


• • n„ ) 

: 


(5.34) 


Hence, if only 4 components of the controller are to be retained they would be 

xj - (^p. f,p. (5.35) 

and if only 2 components of the controller are to be retained, they would be 


- 4 ' 
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® (flpt Hp) 


(5.36) 


The corresponding reduced controller dynamics are deduced from (5.32) in the manner 
described above. 

STEP V-E '. Controller Reduction With Uncertain Parameters (o^. >0) 

Delj^te ^rom^the controller (5.31) those states X- with the smallest component costs 
t/^ = ) computed by 


■'t ■ L ° “ vl ' 


(5.37) 


for j = (''■’■Nr ■’■ 2Nnk), i = 1, ...» Nn« where 


~A 0“1 To 

- • 0 . =■■ 

L p J L P J 


(5.38) 


<2.^ = block diag [(2, ,(3..., o^(2] 

and where Q, is defined by (5.32) 

Pu'^pose of STEP V-B: The component costs (5.37) represent the sum 


V, = V{k.) + I i/((x,)„J 


(5.39) 


where, the total cost is given by (5.15) with 8=0. By use of (5.33) and (5.30) the 
total cost may be written 


R R 

^ = tr = I l^(x., X ) + J l/(x.,x. ). x.| 


(5.40) 


where in (5.37) picks out the first terms in (5.40) due to X,.. The remaining 
terms in (5.40) are due to the states X- and their sensitivities, which are not needed 
for this design "out may be computed from the remaining terms in (5.37). The parti- 
cular range of indices in the U. computation of (7.58) is due to the choice of 
coordinates represented by (5.38). Of course, any other choice of coordinates may 
be chosen with an attendant change in the range of indices jj in (5.37) to identify 
those variables of the augmented state vector X representing controller dynamics and 
their sensitivites , X. and X. . STEP V-B allows sensitivity considerations to influ- 
ence the controller reduction, whereas STEP V-A assumes the parameters are known 
a priori. There is an additional computational burden in STEP V-B and for this step 
assumption A8 must be changed to 

N(. £ 2Np < 2 1. 2N 

due to the fact that (5.37) decomposes into (v+1 ) Liapunov equations each of which 


■ 1-/'/ '1 
\y / /( 


/ /■ -''/‘•vA ' ■ 


is 2NpX2Nj^ [35]. 

Example of STEP V-B: Assume uncertain frequencies Ui , Hence 

= (o)^ .Ug)* =(1. 2) 

= E(uji-oji)^ = Oq = 100 

02 = E{o)2 - ^ 2 )^ = i°0 " 2(100) = 200 

Carrying out the computations (5.37) we find that the controller states arranged in 
their cost order are 




Hj.. n^» ng2» hg2) 


(5.41) 


Hence if only a fourth-order controller is desired we would choose 


"T 


\l ’ *^rl ’ 


(5.42) 


If only a second-order controller is desired we would choose 

Xp = (Pgp Ug^) , < c. 

Controller (5.43) would lead to instability since the rigid body mode has been 
truncated. By comparing (5.34), (5.35) and (5.36), respectively, with (5.41), 
(5.42) and (5.43), the .Influence of parameter sensitivity on the teduced controller 
design can clearly be seen. For example, (5.34) indicates that (n o more 

important than (^gi » *>g-| ) i’' the nominal optimal controller^^ (5. 3^ , whereas (5.41) 
indicates that, from parameter sensitivity consideration, (n , fi ,) is more 
Important. ^ ® 


6.0 CONCLUSIONS 

We have outlined in some detail many of the critical problems associated with 
the control of highly damped flexible structures. The practical problems include: 

• High performance 

* Assembly in space, configuration changes 

• On-line controller software design 

* Lack of test data 

Underlying all of these practical problems is the central problem of modeling 
errors. To Justify the expense of a space structure, the performance requirements 
will necessarily be very severe. On the other hand, the absence of economical 
tests precludes the availability of reliable data before flight. Thus, a gveatei' 
burden of responeibilitij is placed upon analytical methods in the design without 
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the benefit of test data for a sygtem nearly unstable without control. This is 
the triple jeopardy faced by all flexible space structure designs. 

Some precise statements have been verified concerning the performance of such 
structures in the presence of controllers based upon erroneous models. I'he modeling 
errors can always be classified into four categories: (1> parameter errors, (2) 

model order errors, (3) disturbance errors and (4) neglected nonllnearltles. Fin- 
ally, a design algoritlun is offered which has these properties: ‘ 

. 1. Provides damping for a larger numl or of inodes than the optimal attitude 

controller controls. This rate measurement feedback design, with collocated 
rate sensors and actuators, provides control spillover protection for the 
truncated controller 

2. Coordinates the rate feedback design with the attitude control design by 
use of a similar cost function 

3. Provides model reduction and controller reduction decisions which are 
systematically connected to the mathematical statement of the control 
objectives and the disturbance models 

There are many possible versions of the CCA design algorithm. Some versions would 
require more computation with some attendant improvement in performance. For 
example, other choices of coordinates besides "modal" coordinates may provide better 
reduced models in Section 5. However, in the interest of clarity the prtpcedures 
have been explained with some economy of detail. 

Flexible space structure control is indeed a fitting challenge to the best of . 
available control and estimatio.i theories. And, conversely, out of such challenging 
f examples come pointed needs for new theory. '' 
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1.0 INTRODUCTION 


Many proposed large spacecraft have missions which embody two challenging 
requirements: 

(1) Rapid reorientation maneuvers ("slews") 

(2) Rapid vibration arrest and fine pointing upon completion of the slew 

The objective of this discussion is to overview some methods and applications which 
determine optimal maneuver controls. Aside from the usual modeling Issues for large 
space structures, the large-angle maneuver problem is complicated by the inherent 
nonlinearity (l.e., three-dimensional rotational maneuvers are nonlinear even for a 
rigid spacecraft). Of course, additional nonlinear effects arise when even small 
deflections are considered due to kinematic effects such as rotational stiffening. 
Thus, we have the dual curse of nonlinearity and h’gh dimensionality. 

One can approach optimal maneuvers in a fashion analogous to trajectory design 
for launch vehicles. A nominal open-loop maneuv'er could be determined a vviori\ real- 
time feedback controls would then determine small perturbations to the open-loop con- 
trols to account for modeling imperfections and external disturbances. For many 
vehicles and missions, it may prove satisfactory to determine open-loop maneuvers 
which suppress vibration and simply forego feedback controls until the terminal fine- 
pointing phase after completion of each slew. Of course the terminal vibration sup- 
pression and fine-pointing problem is naturally suited to feedback control. In this 
paper, we emphasize determination of the optimal open-loop nonlinear maneuvers 
(slews). The slew controls are determined in such a fashion that vibration of 
selected modes is penalized en route and arrested upon arrival at the ta.get state. 
Some attention is also given to large, slow, linear maneuvers for which the asso- 
ciated slew controls can be cast in feedback form. 

We initially summarize the main aspects of optimal control theory which are 
implicit in the remaining developments. We then discuss the essential ideas involved 
in a class of methods ("continuation" or "homot.opy" methods) which we have found most 
useful in solving the resulting two-point boundary value problems. Several low- 
dimensioned, nonlinear maneuvers of multiple rigid-body configurations using optimal 
momentum transfer are discussed. Sevi-.al linear and nonlinear flexible-body maneuvers 
are then presented and include distributed controls, vibration suppression/arrest, and 
computational issues. Finally, we summarize the status of our present work and the 
key problem areas in which future research appears most urgent. 

2.0 NECESSARY CONDITIONS FOR OPTIMAL MANEUVERS 

We will consider here only finite dimensional systems which can be satisfactorily 
modeled by a system of generally nonlinear ordinary differential equations of the form 


X = f (t, X, u) (2.1) 


where x is the n x i state vector, u is the m x i control vector, and f is 
an n X 1 vector of generally nonlinear functions of all arguments. Suppose, for 
the present discussion, that we are concerned only with the fixeu-end point problem 
for which the boundary conditions 


( 2 . 2 ) 


I 


X (t^) = Xq X <tf) = Xj 

are fully specified. The more general circumstance in which some boundary conditions 
are free or are constrained to lie on a terminal manifold is considered in refer- 
ences 1-6; essentially, some of the 2n conditions (2.2) are replaced by appropriate 
transversality conditions. The number of boundary conditions will remain 2n. Suppose 
we seek an optimal control vector u(t) initiating at x^ and terminating at x^, 

subject to the requirement that a performance functional of the form 

J = I ?'(t, X, u)dt (2.3) 

•'to " ■ 

is minimized. For the present discussion, we restrict u(t) to be continuous with 
two continuous derivatives. The necessary conditions (Pontryagin' s principle form) 
involve the Hamiltonian functional 


H = F + f 


(2.4) 


where A is an n x 1 co-state vector of Lagrange multipliers satisfying the 
adjoint differential equation 



(2.5) 


and the optimal control is determined at each instant by requiring H(t,x,A,u) to be 
minimized with respect to admissible u(t); for twice differentiable f, u, this 
requirement yields the conditions 


8u 



positive definite 


( 2 . 6 ) 


The conditions (2.6) can usually be solved for u(t) to determine the optimal con- 
trol in the functional form 


u = 5 (t,A,x) 


(2.7) 


Finally, substitution of (2.7) into (2.1) and (2.5) allows us to eliminate u and 
obtain the 2nth order system 


X ■ F(t,x,A) x(to> = Xjj x(tf) = X£ 
A = G(t,x,A) 

The above formulation is summarized in figure 1, 


( 2 . 8 ) 








y 
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since the 2n boundary conditions are split, we have a two-point boundary value 
problem (TPBVP) , These differential equations (2.8) are generally nonlinear, often 
exhibit "stiff" character, and generally resist solution except by iterative numeri- 
cal methods. Indeed, a substantial fraction of the optimization literature ("direct" 
methods) have been motivated by the numerical difficulties one often encounters in a 
frontal attack upon equations (2.8). ir. essence, most of the direct methods are 
iterative function space gradient methods which "directly" adjust u(t) and itera- 
tively solve equations (2.1) to minimize (2.3) and satisfy terminal boundary condi- 
tions. Numerical experience suggests that the direct methods efficiency degrades 
rapidly as the problem dimension increases, especially if one does not possess "good" 
starting estimates for u(t); so the issue of which algorithm to apply remains open. 
The recent progress of references 1-6 suggests that improved "indirect" algorithms 
based upon continuation method iterative solutions of equations (2.8) are indeed 
feasible; the methods discussed herein are of this type. The main recent progress is 
associated with greatly decreased reliance upon "sufficiently good starting estimates" 
implicit in the continuation family of algorithms. The essential ideas are indicated 
in the following section. 

3.0 CONTINUATION ItETHODS FOR SOLUTION OF TWO-POINT BOUNDARY VALUE PROBLEfiS 

The basic approach is to embed the problem you wish to solve into a one param- 
eter (a) family containing two prominent members: 

For a = 0, the family reduces continuously to a problem whose solution 
is available 

For u = 1, the family reduces continuously to the problem whose solution 
is sought - ■ 

By sweeping a, one can generate a sequence of neighboring problems (e.g., for 
0 = < nt 2 “3 ^ otfj = 1); by solving each (rxj^) problem in sequence, v;e 

determine "stepping stone" solutions which provide, by extrapolation, starting 
iteratives for the next (uj^^j^) problem. In principle, arbitrarily close starting 
iteratives can be determined. Fxcept for certain singular events (e.g., bifurcation 
points), this approach can very nearly guarantee convergence. Within this approach, 
there exists a variety of ways to construct the f.amily (homotopy chain). The manner 
in which the homotopy family is constructed lias definite irapait upon the ensuing 
algorithms. We have made productive use of three continuation methods which we now 
summarize. 

3.1 Continuation Method 1. Differential Equation Embedding 

The nonlinear system (2.8) is written in the form shown in equation (3.1). This 
form is quite natural for many weakly nonlinear problems; of course, a = 1 is the 
true physical value, and a = 0 results in a linear problem. Provided the linear 
boundary value problem can be solved (which is not a foregone conclusion, especially 
for high degree-of-f reedom systems which possess a wide eigenvalue spectrum, i.e., 
the system is large and stiff), then sweeping a in principle generates a family of 
neighboring problems which can be iterated for ^q('^) to satisfy the final boundary 
condition x(t£) X£, Each Newton Iteration for A^j(a) can be initiated with an 

estimate based upon extrapolating the neighboring a solutions. This method is used 
to solve for optimal maneuvers for several cases in the developments below. 
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The nonlinear state and co-state differential equations (2.8) are rearranged 
and the continuation parameter (a) introduced so tha*" 

X •• X + ® {nonlinear terms} 

(3.1) 

X •= A 2 j^ X + A 22 ^ + ot {nonlinear terms} 

For a = 0 linear problem 

For a = 1 — problem of interest 

Sweep ct to generate a family of neighboring dynamical systems; for each a, iterate 
Xg(a) using Newton's method to satisfy the terminal boundary conditions using neigh- 
boring a solutions for ^o(a) to start each iteration. 

With reference to equations (3.2) through (3.6), we consider the Method 2, 
fixed-point algorithm continuation method of reference 7. This method is called a 
"fixed point algorithm" since it has the property that “ 9 solution of 

the homotopy family at a = 0 for any arbitrary specified This algorithm is 

most useful if one wishes to use an approximate initial guess C for Xq which is 
too far from the solution to permit reliable convergence using Newton's method. 

Method 3 (boundary condition embedding) constructs a homotopy family of 
boundary conditions for those problems which have the property: 

The general nonlinear problem degenerates continuously to a linear 
problem for certain choices on the specif led state bou^dary condi- 
tions. As an example, three-dimensional rigid-body maneuvers are 
nonlinear, but the single-axis maneuver special cases are linear. 

By sweeping a, the specified initial and final boundary conditions are varied from 
the easily solved problem at a = 0 to the problem of Interest at 0=1. 

3.2 Continuation Methods 2 and 3 

Method 2. Fixed-Point Algori thm : 

An algebraic root solving illustration is: 


F(x) = 0 


(3.2) 


Initial guess: x s c 

Construct the homotopy family 


H [zCci),®] = F [z(o)] - (1 - a) F (c) » 0 


(3.3) 


i 



(3.4) 



Note the properties; 


H [z(0),0] 


If we take z(0) = c arbitrary "fixed point" and 


H [z(l),l] H F [z(D] = 0 


(3.5) 


thus z(l) = X by comparison of equations (3.2) and (3.5). 

Method 3. Boundary Condition Embedding; 

When the TPBVP.at hand can be easily solved if the specified boundary conditions 
are given a special set of values, then we can construct the homotcoy family of 
boundary conditions: 


CONTINUATION 


FINAL DESIRED 


"special set OF 

Vector of 


VECTOR OF 


SPECIFIED R.C., 

SPECIFIED 

= a 

SPECIFIED 

+ (1 - a) 

...FOR WHICH THE 

B C 


R C 


TPBVP CAN BE 

» «■ 


13 • Vi* • 

— « 


_ SOLVED 


These considerations are addressed to the important issue of decreasing the 
amount of prior information required to start an iterative solution of equations (2.8). 
For high-dimensional problems with high sei.-itivity to Aq, these methods are most 
important; it may prove Impractical to extract a solutio.i via the classical "shooting 
technique." A separate and fairly well-developed subject matter concerns the cal- 
culation of the successive numerical solutions (i.e., for each a). The most popular 
method is Newton's root solving method in conjunction with Runge-Kutta integrations; 

^ is adjusted based upon linearizing the n boundary condition residuals depen- 
dence upon Aq about the previous iteration. 

Another popular class of procedures is the weighted residual methods. In this 
approach, each state and co-state variable is expanded with undetermined coefficients 
times a finite set of specified functions of time (typically, a set of orthogonal 
polynomials). The coefficients of the series are adjusted by solving a linear system 
of equations obtained from substituting the expansions into the equations of motion, 
linearizing, and evaluating these linear equations at a sufficient number of times 
along the previous trial trajectory. This class of procedures seems particularly 
well-suited to stiff differential equations for which the traditional (e.g., Runge- 
Kutta) numerical methods prove unstable; however, the large linear algebraic cystfms 
encountered are a deterrent to many degrec-of-f reedom applications. 

We mention a third approach, the quasi- linearization method of Miele (the method 
of particular solutions), as developed in references 5, 6, and 8. This approach 
forms successive iterations by: (i) linearizing the departure motion differential 

equations along the previous trial solution of equations (2.8); (ii) determining a 
set of n coefficients to linearly combine n trial integrations of these linear 
equations to produce a new departure motion which satisfies the boundary conditions; 
and (lii) adding the departure motion to the previous iteration to produce the next 
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trial solution. The major disadvantage of this approach is the necessity to store 
each trial trajectory; for many state variables, this storage penalty becomes 
excessive. 

All of the above methods are developed in detail in reference 4 and compared 
using a common set of example problems. Here we note that all of these have been ■ 
used successfully (refs. 4, 6) to calculate optimal nonlinear spacecraft maneuvers; 
if the trial trajectories of equations (2.8) can be integrated accurately using 
Runge-Kutta. methods, then we find the classical Newton iterations used (if necessary) 
in conjunction with a continuation method (to provide accurate starting iterations) 
to be preferred. However, for very stiff differential equations, one must resort to 
other methods such as those discussed above, or Indeed, abandon the indirect approach 
altogether.. It is apparent, however, that a large fraction of the optimal nonlinear 
maneuver problems of practical interest can be solve d by the indirect methods 
illustrated herein and documented in the references. 

4.0 OPTIMAL LARGE-ANGLE MANEUVERS OF MULTIPLE REACTION WHEEL SPACECRAFT 

With reference to figure 2, we consider the dynamics of a generally asymmetric 
vehicle B having three identical rotors aligned with the principal axes. The vehicle 
is maneuvered via internal motor torques [uj^(t), U 2 (c) , u-j(t)3 applied to the wheels, 

[-Uj^(t), , -u^(t)] applied to assumed rigid body B, and external torques 

[Lj^(t), L2(t), L^(t)3 from unspecified actuators. The equations of motion are 
equations (4.1) - (4.3). 


Kit^mat ics : 

O 

(wj^Bq — W2B3 u)267)/2 
6^ = (■^^1^2 u)2Bj^ + 0)380) /2 



(4.1) 


Spacecraft Dynamics : 


-1 



(I3 — 12)^20)3 + h30)2 


h20)3 





ac 












(4.2) 
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Wheel Dynamics ; 




1 - 1,2,3 


(4.3) 


4.1 Optimal Maneuver Necessary Conditions for the Four-Body Configuration 
Performance Indices ; 

■’l ■ i '■2 '■D "2("1 + “2 + "3)] ■ X'' ‘1 

^ /o'* ^3)-" "2(“' * “2 “3)] ^ ■1' 

■’2 ■ i /o'* ^3) “2(=? * “2 * “3)] ■“ ■ /o'' ^ ■*' 


F4 <it 


(4.4) 


(4.5) 


(4.6) 


(4.7) 


Hamiltonian: 


T • T • T • 

H = Y e + X U) + a h + F, 


(4.8) 


where the and X^^ are Euler parameter and angular velocity co-states, and 

the Oj are wheel relative momenta co-states. 

Pontryagin Necessary Conditions ; 


Minimize H with / 
respect to Uj^(t) S 


^-0 0- 1,2,3 


3H 


for of equation (4.4) (4.9) 


- 0 i - 1,2,3 
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Co-State . 

differential 

equations 


r 

dt 

3H 
“ 3B^ 

i - 

0,1,2. 3 

(4.10) 

dX^ 

dt 

3H 

" 3^1 

i - 

1,2,3 

(4.11) 

^ dt “ 

3H 

■ 3h^ 

1 = 

1.2.3 

(4.12) 


(+ boundary conditions) 

Notice we have a once redundant ten-element state vector 




^2 i “l ^2 “3 J ^1 ^2 ^3 


ii: 


w 


where 


(4.13) 


6 is the Euler parameter vector which describes the inertial orientation of B 
w is the inertial angular velocity of B 

h is the momentum of the reaction wheels relative to B (h “ 0, for locked 
wheels) - o 


The Euler parameters satisfy the constraint 

b'^ 6 ° 1 


(4.14) 


Their use makes the formulation universal, as opposed to the transcendental differen- 
tial equations with a singularity which one obtains for any choice of Euler angles. 
However, one often requires a transformation from Euler angles to Euler parameters; 
for the 1-2-3 set of Euler angles (Sj^, 02 j 63 ) » the transformation is 


Bq - c(01/2) c(02/2) c(03/2) - s(01/2) s(02/2) s(03/2) 

- s(01/2) c(02/l) c(03/2) + c(0l/2) s(02/2) s(03/2) 

&2 - c(01/2) s(02/2) c(03/2) - s(01/2) c(92/2) s(B3/2) 

B 3 - c(01/2) c(02/2) s(03/2) + s(01/2) s(02/2) c{03/2) 

c( ) = cos ( ), s( ) H sin ( ) 


(4.15) 


f 












// 
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These results and similar developments for all twelve classical Euler angle sets 
are given in reference 4. 

We. must now consider the necessary conditions [in addition to satisfying, equa- . 
tions (4.1) - (4.3)3 optimal maneuvers of the configuration of figure 2. First, 
one must define how to measure optimality. Equations (4.4) - (4.7) define four' 
integral measures of optimality (from refs. 4, 5, 6), The first of these can be 
Interpreted as a positive measure of rlie control effort, while the remaining three 
penalize control derivative variations. As one might expect, the derivative 
penalties lead to smoother optimal torques with lower frequency content; such con- 
trols are attractive wtien generalized to flexible vehicles because they tend to excite 
only lover modes (which arc usually better modeled). 

.Notice the ten co-state equations couple together with the ten state ecoatious 
(equations (4.1) - (4.3)] to yield a 20th-order system. Typically, the ten initial 
co-ctates must be determined Iteratively to satisfy ten constraints at the final 
t Ime . 


In order to discuss solution techniques, it is useful to consider particular 
maneuvers. Four cases are now discussed. Table 1 summarizes some macroscopic 
Information concerning the four cases. Case 1 is the simplest; the solution of 
this particular case is treated in detail in reference 2. Since the wl\eels are 
locked, the system becomes a single rigid body with seven state variables 
xT “ [Sg S], ^2 ^3 ‘*^1 '^2 '^ 3 ]- reader can readily verify that any set of "pure 

spin" boundary conditions causes all of the gyroscopic nonlinearities of equa- 
tions (4.1), 4.2) and (4.10), (4.11) to vanish and results in a linear system which 
is readily solved analytically; reference 2 gives the details. Thus, "boundary con- 
dition embedding" is motivated as an attractive way to solve the TPBVP by 
continuation. " o ° 

0 , 

4.2 Four Maneuver Cases 

*^or this case, table 2 gives the corresponding Euler angle boundary conditions 
and equation (4.16) gives the boundary condition embedding (homocopy chain) family 
of boundary conditions: 


®10 

(ct)" 


0 


0 

®20 

(a) 


0 


0 

®30 

(a) 


0 


0 

‘"10 

(a) 


0.01 


0.01 

“20 

(a) 


0.005 


0 

“30 

(a) 


0.001 


0 

— 

— 

= a 

— 

+ (1 - a) 

— 

®lf 

(a) 


511/2 


511/2 

®2f 

(a) 


TT/3 


0 

CD 

J'ti 

(a) 


ir/4 


0 

“if 

(a) 


0 


0 

“2f 

(ct) 


0 


0 

_“3f 

(a) 


0 


0 


(4.16) 


We found 02 , ct^} = {O, 0.0001, 0.25, 0.5, 0.75, l} resulted in reliable 

convergence for large variations in boundary conditions. 


We note in passing that the a = 0 case is a pure spin single-axis rotation 
about the b, axis; this case can be solved analytically (ref. 2). We performed 
other continuations corresponding to rotations about b 2 and b^ and obtained the 

same final solution. In this case, it was found that introducing the a-sequence 
{ 02 ^... ag} = { 0 , 0.001, 0.25, 0.50, 0.75, 1} and using Newton's method to iterate the 

initial co-states (Xi(to) X 2 (tQ) X 3 <to) Yi(tQ) 'Y 2 (to) Y 3 (tQ)} = converged 

efficiently. Specifically, the iteration for initiated using the fol- 

lowing starting procedure: 


For tt 2 = 0.001, Aq(® 2) “ ^o(°^l^» calculated from analytic solution 

= 0.25, ^ 0 ( 0 ^ 3 ) extrapolated linearly from the converged 

“ 4 * O 3 , otg, extrapolated quadratically from the converged 

^o <“!-!> • ^o<“i-2>» ^oK-3> 

2 

As is evident from table 3, the convergence tolerance was held loose (10~ maxi- 
mum relative error of final boundary conditions) until the final (Og = 1) continua- 
tion. This iteration process proved rather efficient as can be judged by the small 
number of Newton iterations. The history of state and control variables on the con- 
verged optimal maneuver is graphed in figure 3. Thus the final maneuver, which does 
not resemble a single-axis maneuver, is easily obtained via a continuation process 
initiated with the single-axis special case. 
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We now consider Case 2 (table 1). For this case, only the wheel along b£ is 
free -'-o rotate. The initial state is a “flat spin”; the body rotates uniformly 
about bg (the largest axis of inertia); it is desired to transfer all angular 
moaientum to the wheel via a judicious motor torque u(t). This problem has received 
considerable attention historically (e.g., refs. 9, 10) for^the constant torque 
(u(t) = constant) and constant relative momentum transfer (h(t) = constant) case. 

In figure 4 (the right half of the figure), we display the constant torque 
(u(t) = 0.005 N*m) maneuver published by Barba and Aubrun (ref , 9) . The general 
features of the maneuver are the following; 

(1) The maneuver is nonlinear; a separatrix occurs around 1000 sec prior to 
which the body b^ axis cones about the inertlally fixed momentum vector 
and after which the wheel axis (b 2 ) cones about the fixed momentum vector 

(li) Without damping or another torque mechanism, the momentum transfer 

maneuver cannot be perfectly realized; i.e., the body cones at a nonzero 
final mutation angle 0 between b 2 and the Inertially fixed angular 
momentum vector 

We address the issue as to whether a nonconstant torque would significantly improve 
this flat-spin recovery maneuver. As is evident in the left half of figure 4, the 
optimal maneuver is qualitatively similar to the Barba/Aubrun maneuver, but the final 
mutation angle is decreased by 50 percent. The optimal torque history fluctuates 
±10 percent about 0.05 N’m during the first 1000 seconds as is evident in figure 4. 

The flat-spin recovery, while Interesting, does not represent a very general 
class of maneuvers, since only one final rest orientation is stable (and it cannot be 
reached exactly without bringing other control or energy dissipation elements into 
action) . Cases 3 and 4 are provided to demonstrate maneuvers which bring a generally 
rotating body to a state of rest in a general final orientation. These maneuvers 
also demonstrate the practicality of calculating highly nonlinear maneuvers with up 
to six control variables and of solving up to a 20th-order system of simultaneous 
nonlinear equations. 

Case 3 is a generalized optimal momentum transfer using all three wheels simul- 
taneously. In essence, the objective is to transfer all of the angular momentum intc 
the three wheels and leave the body at a prescribed inertial pointing. Figure 5 
shows the state and control variable history. Notice, since this is an internal 
torque maneuver, the total angular momentum is fixed; therefore, the final wheel spec' 
can be calculated from the initially calculated angular momentum. 

Case 4 is a maneuver involving both internal and external torques. The 
criterion 
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Is minimized. The terminal boundary conditions- require that the body be at rest and 
all wheel speeds be reduced to zero. Thus this maneuver has the simultaneous objec- 
tive of dumping the system momentum. As is evident in figure 6, the maneuver is 
quite smooch and the terminal boundary conditions are satisfied. Cases 3 and A con- 
trast sharply with Case 2; the nutational motion associated with a one-wheel configu- 
ration is not evident due to (1) the lower wheel speeds, and (ii) the ability to 
generate optimal torques in any desired direction. 

The foregoing results clearly suggest that we can, at least for a wide spectrum 
of problems, routinely determine optimal nonlinear maneuvers for spacecraft having a 
moderate number of degrees of freedom. Of course, nonlinearity is only one compli- 
cating issue for large flexible spacecraft; to what extent -can we generalize these 
applications for vehicles where flexibility is a paramount concern? Specifically, 
can we apply similar methods to maneuver flexible vehicles, penalizing vibratory 
motion en route and arresting certain modes upon arrival? As is evident in the 
following developments, we have made some significant progress toward answering 
these questions in the affirmative. 

5.0 SLEWING MANEUVERS FOR. FLEXIBLE SPACECRAFT 

The optimal control problem for maneuvering a vehicle of the shape shown in 
figure 7 is complicated by the presence of the attached flexible appendages. Mathe- 
matically, the presence of flexibility in the plant description Increases the number 
of degrees of freedom required to model the motion of the vehicle (refs. 1, 3, A, 
and 11). As a result, in order to successfully maneuver the vehicle, the control 
designer must have an accurate plant description which accounts for both rigid and 
flexible effects. In addition, when high angular rates are achieved during the 
maneuvers, great care must be exercised in modelling the vehicle and formulating 
systematic algorithms for solving the resulting nonlinear optimal control problems. 
Hov;ever, the Issues of model truncation and control spillover effects are not dealt 
with in this paper. 

Section 6 develops an open-loop optimal control strategy for linear large-angle 
single-axis maneuver of a flexible spacecraft possessing a distributed control 
system. Section 7 presents a closed-loop formulation for maneuvering flexible 
vehicles. Section 8 develops a differential equation embedding continuation method 
for solving nonlinear open-loop formulations of the flexible spacecraft equations of 
motion. Example maneuvers are given in section 9. 

,6.0 LINEAR LARGE-ANGLE MANEUVERS OF FLEXIBLE SPACECRAFT USING DISTRIBUTED CONTROL 
6.1 Equations of Motion 

The linear time— invariant form of the equation of motion for the vehicle 
depicted in figure 7 Is given by (see refs. A and 11) 


MS + K5 = Pu 


( 6 . 1 ) 
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where 



P = 


1 

0 





[ 1 , 1 1 ] 


(N xl) 
ac ' 


[F] 


ij 




- r 


1 . 


n 



N 


ac 


where n is the number of flexible modes in the model, N__ is the number of 
appendage controls, x^ denotes the point of application of the jth appendage con- 
trol, 6 is the rigid-body angle, n is the vector of modal amplitudes, Uj^ is 
the rigid-body torque, U 2 is the vector of appendage control torques, r is the 
radius of the rigid hub, the integral definitions of M 0 jj, are given 

in references 4 and 11, M is positive definite, and K is positive semidef inite. 
In equation (6,1) it has been assumed that quadratic terms in the modal amplitudes, 
amplitude rates, and angular velocity are vanishingly small and thus deleted. 

For the optimal control problem, equation (6.1) is transformed to modal space 
and cast in the first-order form given by 


s = As + Bu 


where 


( 6 . 2 ) 


A = 



I 

0 




e'^me = I 


T 

E KE 


A = diag(0, A^, 


..., A^) 


where E is the normalized eigenvector matrix for M and K, and A is the diagonal, 
matr.lx. containing the eigenvalues of M and K. 
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6.2 Optimal Control Formulation 

We seek an optimal solution of equation (6.2) which first satisfies the 
prescribed terminal boundary conditions given by 


r T. “|T 

r. , 

•T, "|T 

(6.3) 

So “ n ^o>J 

So = L9^*=o> 


Sf = n’^<’'f^] ^ 

y “ (j(tf) 


(6.4) 


and second, minimizes the performance index given by 


J = ~ f ^ fu’’^ W u + s'*’ W si dt 
2 Jt^ L- uu - - ss-J 


(6.5) 


The selection of the performance index of equation (6.5) is arbitrary although made 
for convenience. The elements of the state weighting matrix Wgg are selected to 
make the state penalty term in J proportional to the kinetic and potential energy 
of the vehicle. Since it is desirable to suppress the elastic deformations and 
deformation rates at the end of the slew, we also impose the constraint 
n(tj) = n(t£) “0 in equation (6.4). 

The control weight matrix is chosen to be diagonal with smaller weighting on the 
rigid-body control torque than on the appendage control torque. This is because the 
elastic appendage controllers are meant to serve as vibration suppressors while the 
rigid hub controller is meant to execute most of the slewing maneuver. The penalty 
on the rigid-body angle is also chosen to be small because the structure is required 
to undergo a large-angle rotation. The rigid-body angle penalty Is not set to zero 
since this choice leads to numerical problems in computing the state transition 
matrix. However, by adjusting the rigid-body angle penalty to be two to three orders 
of magnitude smaller than the penalties on the modal amplitudes, amplitude rates, and 
angular rates, the potential numerical problems can be overcome. 

Using Pontryagin's principle (ref. 4), the necessary conditions defining the 
optimal maneuver are given by equations (6.6), (6.7), and (6.8): 

State: 


s “ As - BW 
- - uu - 


( 6 . 6 ) 


Co-state: 


A = -W s 
- ss- 


T 

- A^A 


(6.7) 


Control: 


— 1 T 

u « -W B A 
- uu - 


( 6 . 8 ) 
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The solution for equations (6.6) and (6.7) Is given by 

S(tf) 


where e 


x(tp 


«(tf-to) 


n(tf-t^) I 


M^o) 


ft 


-1 T 

-BW 

uu 


-A* 


(6.9) 


is the exponential matrix. Since s(tp) and s(t^) are known and 
X(t ) and A(tf) are unknown in equation (6.9), we need to solve for A(t ) to 

obtain the complete solution. Upon setting ()> * e i o the partitioned matrix, 
solution for s(tj) is given by 


s(tp = <l>sgS(t^) + 


( 6 . 10 ) 


Solving equation (6.10) for the initial co-states yields 


i(to> = [<!>sX^'^ts(tf) - ^^^^(t^)} 


( 6 . 11 ) 


The optimal control time histories follow upon integrating equations (6.6) and (6.7) 
subject to s(t|j) given by equation (6.3) and X(tp) given by equation (6.11). 
Example maneuvers are given in section 9. 

7.0 LINEAR LARGE-ANGLE MANEUVERS OF FLEXIBLE SPACECRAFT USING FEEDBACK CONTROL 

7.1 Equations of Motion 

The equations of motion given by equations (6.1) and (6,2) are used for the 
plant description of the feedback control problem of this section. 

7.2 Optimal Control Formulation 

The optimal feedback control for slewing a flexible spacecraft differs funda- 
mentally from the open-loop problem of section 6.2 in two ways. First, the instan- 
taneous values for the state influence the applied control. Second, the boundary 
conditions for the maneuver are specified initially and are free at the final time. 

Breakwell (ref. 12) has recently presented two approaches for feedback control 
of a flexible spacecraft. First, he presented a constant feedback gain approach. 
Second, he developed a time-to-go formulation with time varying gains. The perfor- 
mance index he selected is of the form 


-f 


SSf 




dt 


(7.1) 
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where S Is the terminal state weighting matrix, Q Is the state weighting matrix, 
and R Is the control weighting matrix. Breakwell (ref. 13) has also developed a 
dlstrlbuced transfer function approach for the feedback control problem; however, 
his distributed transfer function approach Is not presented in this paper. 

Breakwell's constant gain solution Is given by equation (7.2) and his. time varying 
solution is given by (7.3): 

tj » «■; u = -r“^b'^P^ (7.2) 

« oo; u - R"V[s"^ (|) 22 (t) - <f>i2^'^^~^[5ll^'^^ " '*’21^^^ 5^*=^ 

where Is the solution to the algebraic matrix Rlccatl equation and T Is the 

time to go. Example maneuvers are given in section 9. 

8.0 NONLINEAR LARGE-ANGLE MANEUVERS OF FLEXIBLE SPACECRAFT USING DISTRIBUTED CONTROL 

8.1 Equations of Motion 

When high angular rates are achieved during slewing maneuvers, the equation of 
motion given by equation (6.1) must be modified in two ways. First, quadratic terms 
in the angular velocity are retained in the equation of motion. Second, an arc length 
correction is added to the kinetic and potential energy integrals in order to 
properly account for second-order effects in the mass and stiffness distributions. 

As a result, the nonlinear equation of motion given by equation (8.1) accounts for 
the so-called centrifugal stiffening effect of the elastic appendages: 


MC + = ? u 


where 


( 8 . 1 ) 



The Integral definitions of M„ , M , M*, and K ^ can be found in references 4 

and 11. 
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Equation (8.1) is transformed Into modal space and cast In the first-order form 
for the optimal control problem using the coordinate transformation defined In 
equation- (6.2) : 


8 “ A(s,a) s + B u 


where 


A(s,a) = 


A2j^(s,a) 0 


A2i(s,a) “ -A - L 


( 8 . 2 ) 



and a Is the continuation parameter which has been introduced to aid In the solu- 
tion of the nonlinear problem. 

8.2 Optimal Control Formulation 

We seek an optimal solution of equation (8.2) which satisfies the prescribed 
terminal boundary conditions given in equations (8.3) and (8.4) and minimizes the 
performance index of equation (8.5). 



(8.3) 

(8.4) 

(8.5) 


As part of the prescribed terminal boundary conditions in equation (8.4) we also 
Impose the constraint that n(tj) = n(tj) = 0, 

Since equation (8.2) is nonlinear, a closed-form solution for the initial 
co-states is not possible. However, the artificially introduced continuation param- 
eter in equation (8.2) can be varied from 0 to 1 in a manner which embeds the non- 
linear problem in a family of problems which includes the linear problem of section 6 
(l.e. a “ 0) as a special case. Thus the solution for a = 0 can be used to 
establish a sequence of problems which converges to the nonlinear problem of interest 
(namely, when a * 1) . 




The necessary conditions defining the nonlinear optimal solution are given by 
equations (8.6), (8.7), and (8.8). 


State: 


Co-state: 


8 - A(s,a)s - BW"Vx 
- uu - 


( 8 . 6 ) 


X “ -W s - C(s,a)X 

““ S S“ “ 


where 


C(s,a) 


(8.7) 


T 

A2^(s,a) 


T T T 

I -2a(e^S2)esJ^L 


Control: 


-1 T 

u = -W B X 
— uu - 


( 8 . 8 ) 


The partial derivatives of equations (8.6) and (8.7) required In the nonlinear state 
transition matrix are given in references 4 and 11. 

8.3 The Continuation Method Using Differential Equation Embedding 

For the solution of the nonlinear optimal control problem, we introduce the 
sequence of continuation parameters {O = < ... = l), where p Is either 

preset or determined during the solution process. The operator equation defining the 
optimal control solution is given by 



s — A(s,a.)s + BW ^b"^! 
- 1 - uu - 

1 

( > 
0 


X + + C(s,a^)X 

> = { 

l- 


2f (desired) " V '^i<>^o> »^"i)j 


bJ 


(8.9) 


The continuation method of this section seeks the solution equa- 
tion (8.10) using extrapolated estimates of based on back a-values. 




( 8 . 10 ) 


The procedure is an iterative method since F is nonlinear; however, we have found 
that only two to four Intermediate a-values are typically required to solve the 
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.‘t -I . 


nonlinear problems of this section. In the nonlinear examples of section 9 the 
iterative solution for each a-value required at most two or three Iterations to 
converge. 

The differential correction strategy is to seek the correction vector AX 
subject to the terminal constraint given by 

Sp/j 4 - s(x + AX, o = d X =» X(t ) (8.11; 

-f (desired) -'-o — ’ f' - -o - o' 

for a specific a-value. Upon linearizing equation (8.11) we obtain equation (8.12), 
where s denotes the numerically integrated solution of equation (8.6) using the 
approximate initial co-state X^ 


-f (desired) - 3X 

~o 


AX = 0 


The equation defining the solution for ^ is then given by 


ax . 

-O tf 




— f -f (desired) - 


The solution for M in equation (8,13) is easily obtained using Gaussian 
elimination. 

For each a-value in equation (8.10), equation (8.13) is iteratively solved 
until the norm of is less than some small value, that is < £. Then a 

is incrementally increased, and the solution for the intermediate a-value is 
obtained from equation (8,13). The process continues until a = 1 and the non- 
linear problem of interest has been solved. The numerical algorithm is summarized 
in figure 8. 

9.0 EXAMPLE MANEUVERS 

Table 4 summarizes the boundary conditions for the example maneuvers of this 
section. The structural parameters for the model used in Cases 1, 2, 3, 6, and 7 
are given in references 4 and 11, and those for Cases 4 and 5 are shown in refer- 
ence 12. 


Cases 1 and 2 (figs. 9 and 10) demonstrate that by adding additional controls 
on the structure the system performance can be improved. In this particular case, 
the first mode peak amplitude is decreased by 9 percent and the rigid-body peak con- 
trol torque is decreased by 4 percent . 
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Case 3 (fig. 11) presents a ten-mode case using distributed control and 
demonstrates that all ten modal amplitudes and amplitude rates (not shown) have 
satisfied the prescribed boundary conditions. 

Case 4 (fig. 12) presents the constant feedback gain solution of reference 12 
and indicates at least for the one-mode case that the vehicle's controlled perfor- 
mance is not unreasonable. 

Case 5 (fig. 13) presents the time-varying feedback gain solution of refer- 
ence 12. We see Chat the boundary conditions are satisfied and that for large 
the control time history is very similar to the open-loop solution. 

Case 6 (fig. 14) presents a somewhat counter-intuitive spin-reversal maneuver 
where the vehicle backs up before moving forward. It can be shown (refs. 4 and 11) 
that the spin-reversal phenomenon is the result of fixing either the maneuver time 
or the final maneuver angle. 

Case 7 (fig. 15) presents a stressing spin-r-revei-sal maneuver where the peak 
modal amplitudes for the nonlinear maneuver differ from the linear solution by about 
17 percent to 41 percent. In particular, the first mode had small amplitudes while 
the higher modes had larger amplitudes when compared with the linear solution. The 
shape of the modal amplitude responses is also slightly different from the linear 
solution. 
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TABLE 1. - FOUR LARGE-ANGLE MANEUVERS OF THE VEHICLE IN FIGURE 2 


CASE I OCSCftIPTION 


rCRfORNMCC 

IHDCX 


SYSTEM 
PARAMETERS 


iMtTlAL 

COHOITIOHS 


DESIRED 
MHAl CONDITIONS 


1.0NTINMTI0N 

method 


mjHERICAL METHOD 

TO SOIVC ODEt 


CENTRAL MANEUVER 
OF A SINGLE 
RIGID BODY 


il\ ♦ l| ♦ 


STATE: B.w 
COfFTROL: L 




*86.215 




o' 

h 

= 

BS.07 


w 

- 

0 

h 


L13.565 




0 


WHEELS LOCKED 




’.01 r/s 
|.DOS 
LOOl 
0 S 


|?3^V 
w(t,) . Q 
tf - 100 S 


5«/2 

■/3 

«/< 


BOUNOAAT 

COt^ITION 

lAEDOlNG 


HEUTOCI DIffEREHTIAL 

corrections of suc- 
cessive 
ruwc-kutta 

SOLUTIONS 

ALSO SOLVED USING 

VE1GKTE0 residual 
METHOD (REF. S) 

"METHOD OF PARTICU* 
lAR SOLUTIONS* 
[Tin, S ANO €) 


optimal 

FLAT SPIN 
RECOVERY* 

USING ONLY THE 
PITCH WHEEL 

STATE: hy 

CONTROL: u. 




♦ [hjit,) - h]* 


SAME AS ABOVE 
EXCEPT -2* 

WHEEL IS 
NOT LOCKED; , 
■ 3 .05 kg wr 

H - I3 «j(V 

• 20 «fS 


w(t„) 


.176 r/t\ 

w - ° 

L • 0 S 


w(tf) • 0 

" H 

If • 4000 S 
IKPOSEO 

TRANSVtRSATlim 

CONDITIONS: 

A(t^) -0 

g(ty) ■ 0 


FIXED 

POINT 

ALGORITHM 


NEWTON DIfFEPENTIAL 
COAfiECriONS OF SUC- 
CESStVC 
RUNGC*KUTTA 
SOLUTIONS 

ALSO SOLVED USING 

*ICTHOD$ OF PARTI. 
CULAR SOLUTIONS* 
(REFS. S AND «) 


OPTIMAL INTERNAL 
TORQUC MANEUVER 
USING 3 WHEELS 




STATE : 

CONTRX: u 



■> 


86.215 

♦ u| ♦ l!|) dt 


■ 

8S.07 


,*3 


113.5G5 


kg« 


Jf ■ 3 • .OS kg ** 
H • INITIAL ANG. 






.6427BD 

.442275 

.442275 

.442275 

.01 r/» 

.005 

.001 




•II 


u(t^) • 0 
h(t^) - H 


KIELE'S -HETHOO OF 
PARTICULAR SOLUTIONS*, 
QUASH IHEARIZATION 
KETHOU» AS 
OtVELOPED IN 
REFS S AND 6 


[optimal 

IMfR:iAL/tXTLO*.Al 

■TORQUE 

Kueuver 


STAH : B.u.h 
CONTROL: u.l 


r (l' 


■4' (>■? ‘ A 


uj) dt 



10 kg m2 

■ 

6.33333 


9.U66A 

• 

J • 5 hg 


RlO 


l(*o) 

hit,) 


.01 r/s 
.005 
.001 

5 KHS 

2 

3 


f.43047 

^.'.70106 

.09230 

1.56098 

u(t^) • fl 

h{lf) - fi 


MIELE'S “METHOD OF 
PARTICULAR SOLUTIONS* 
O'JASlLP.LARIZAtlON 

hlthoo, as 

OevEIOPEO IK 
REFS S AND 6 




TABLE 2,- 

RIGID BODY MANEUVER 

(CASE 1) 


dl. 

Igi I 3 ) 

(1, 0.833333, 0.916667) xlO® kg 





DESIRED FINAL STATE 



INITIAL 

STATE (tjj = 0 s) 

(t^ = 100 s) 



- e,(t^) = 

®io = o 

®l(^f) “ ” &it /2 

1-2-3 Euler 
Angles 

< 

w - 

® 20 =® 

egCtf) = 62 f = tr/3 



, ®3<^o) = 

®30 = '^ 

W ®3f ’ 



- 0 )j(t^) = 

“10 ~ 

0 

II 

3 

11 

3 

Orthogonal 
Ang. Vel. 
Components 

< 


ui^Q = 0.005 

W 2 (t^) = (Dgf = 0 



^ - 3 (^ 0 ) = 

W 30 = 0.001 

W 3 (tf) = ( 1 ) 3 ^ = 0 



TABLE 3.- 

BOUNDARY 

CONDITION 
[From ref . 

CONTINUATION FOR 
A] 

CASE 1 




0 

0.001 

0.25 

0.50 

0.75 

1.00 


Tolcra^cc■ — 

... 

10-* 

10-2 

10-2 

10-2 

10 


No. of iicralions — 


1 

4 

8 

6 

4 


Iniiia! 

•lOo) 

0 

0 

0 

0 

0 

0 

sute 


0 

0 

0 

0 

0 

0 



0 

0 

0 

0 

0 

0 


U,Ug) 

0.01 rad/$ 

0.0100 

0.0100 

0.0100 

0.0100 

0.01 


Uj(fg) 

0 

0.00005 

0.00125 

0.00250 

0.00375 

0.005 


<->j ('o ) 

0 

0.00001 

0.00025 

0.00050 

0.00075 

0.001 


Final 


5r/2 

5i/2 

5t/2 

5t/2 

5ir/2 

7.85397 

Slate 

0 

0.01048 

0.260829 

0.52249 

0.78471 

(-5i/2) 

1.04731 



0 

0.07857 

0.19620 

0.39248 

0.58909 

(-X/3) 

0.78546 


«f <'/> 

0 

0 

0 

0 

0 

<-»/4) 

0 


»2Uf) 

0 

0 

0 

0 

0 

0 


"3 C/* 

0 

0 

0 

0 

0 

0 


Co.nverged 

(fff) 

-4.3124 

-4.31237 

- 4.22029 

-4.01835 

-3.81454 

- 3.64945 

co-siate 


0 

-0.01576 

-0.37601 

-0.67991 

-0.91423 

- 1.10589 

(xlO*) 


0 

-0.03034 

-0.73721 

- 1.34988 

-1.80119 

-2.13870 


T«(7«) 

0 

. 0 

0 

0 

0 

0 


># ('•) 

-0.17649 

-0.17649 

-0.17685 

-0.17546 

-0.17082 

-0.16634 



0 

-0.00063 

-0.01501 

-0.02674 

-0.03496 

-0.04117 



0 

-0.00039 

-0.01101 

- 0.02672 

-0.04411 

- 0.05838 

*MaRifnum teUiivc error in *#’s (noimali/cd by 0.0011 
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TABLE 4.- DESCRIPTION OF TEST CASE MANEUVERS 


Qualitative Description 

No. of 
vndes n 

Rest'tO'Rest Hanuever 
60 sec 

2 

Rest-to-Rest Hanuever 
> 60 sec 

2 

Rest-to-Rest Hanuever 
■ 5 sec 

10 

Rest-to-Rest Hanuever 
t^ * 5 sec 

1 

Rest-to-Rest Hanuever 
Time Varying Feedback 
« 5 sec 

1 

Spin-Up Maneuver 
t^ ■ 60 sec 

3 

Nonlinear Spin Reversal 
t^ » 60 sec 

3 


red/s 


if, rad/s 


I = diag(10"’. 1 
T ■■ diag(10"^, 1 


(2(ntl) : 2(ntl)) 
























PHYSICAL PROBLEM 


PONTRYAGIH’S PRIHCIPlE 


TWO-POINT-BOUNDARY 
VALUE-PROBLEM (TPBVP) 


SYSTEM STATE ODE'S 
X = f.{t,x,u) 

FIND OPTIMAL CONTROL u{t) 
TO MINIMIZE 

0 = ♦(x.) + f J (t,x,u)dt 

SUBJECT TO (FOR EXAMPLE) 
THE BOUNDARY CONDITIONS 

= x^ 

x(t^) = Xf 


HAMILTONIAN 
H = / * f 

CO-STATE EQUATIONS 



MINIMIZE H w.r.t. 
ADMISSIBLE u(t) AT 
EACH INSTANT - THIS 
YIELDS 


FIND ^ SO THAT THE 
SOLUTIONS OF 
STATE A - F(t,x.X) 

CO-STATE X = G(t,x,X) 
INITIATING AT 

■ io 

i(to) = Xq 

TERMINATE AT 
t ^ 


Figure 1.- The optimal maneuver necessary conditions and the associated TPBVP. 



ASYMMETRIC 
MAIN BODY B 


veloelt# ^ 


/ Orbit 
/ Norful 

Wheel Momenta: h^ = J U-; 

i = 1,2.3 

B's Angular Velocity: (10^,0)2,102) 

Euler Parameters 
orienting B: (Bq.Bj^, 32*23) 

Internal (Motor) 

Torques: (u,,ii2.U2) 

External Torques: ;L^, 12,12) 

Figure 2.- Dynamics of a four-body configuration with internal and external torques 
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Optimal maneuver of Vadali & Junkins [ref. 5] 


Maneuver of Barba & Aubrun [ref. 9] 

























Figure 7.- Modal structure. 


ClVeK; s(0).s(tf) 


APPROXIMATE »(0) FROM £Q. (6.11) 


INCREMENT a 


SOLVE THE DIFFERENTIAL EQUATION OF EQ. (8.6) and (8.7) TO DETERMINE 
S ■ s(A(0),tj) 

AND DETERMINE THE PARTIAL DERIVATIVES OF EQ. (8.12) 


CALCULATE THE RESIDUAL VECTOR AS, 


VIA QUASI-NEHTI 
CALCULATE 


PIf 


HODS 


IF SMALL AND a - 1 END 

IF SMALL AND o P 1— A 

IF LARGE. CONTINUE 

I 

calculate ^ FRtW EQ. (8.13) 

I 

APPLY CORRECTIONS 

X(0) - X(0) 4 ^ 

1 

B 

Figure 8,- Differential correction algorithm for the differential embedding method. 
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